автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Математическое моделирование процессов высокотемпературной релаксации в электронно-атомных системах

кандидата физико-математических наук
Прибиш, Ян
город
Дубна
год
2005
специальность ВАК РФ
05.13.18
Диссертация по информатике, вычислительной технике и управлению на тему «Математическое моделирование процессов высокотемпературной релаксации в электронно-атомных системах»

Автореферат диссертации по теме "Математическое моделирование процессов высокотемпературной релаксации в электронно-атомных системах"

(

ОБЪЕДИНЕННЫЙ ИНСТИТУТ ЯДЕРНЫХ ИССЛЕДОВАНИЙ

11-2005-142

4 На правах рукописи

УДК 519.6, 517.9 533.9, 539.18

ПРИБИШ Ян

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ВЫСОКОТЕМПЕРАТУРНОЙ РЕЛАКСАЦИИ В ЭЛЕКТРОННО-АТОМНЫХ СИСТЕМАХ

Специальность: 05.13.18 — математическое моделирование, численные методы и комплексы программ

Автореферат диссертации на соискание ученой степени кандидата физико-математических наук

Дубна 2005

Работа выполнена в Лаборатории информационных технологий Объединенного института ядерных исследований

Научные руководители:

кандидат физико-математических наук Э.А. Айрян

кандидат физико-математических наук Б.Ф. Костенко

Официальные оппоненты:

доктор физико-математических наук, профессор Л.А. Севастьянов

кандидат физико-математических наук А.Е. Волков

Ведущая организация:

Тверской государственный университет

Защита состоится "18 " г. в " "на заседании диссерта-

ционного совета Д. 720.001.04 в Объединенном институте ядерных исследований (Лаборатория информационных технологий), г. Дубна Московской области.

С диссертацией можно ознакомиться в научной библиотеке Объединенного института ядерных исследований

Автореферат разослан 47 " окГл %а*2005 г.

Ученый секретарь диссертационного совета

кандидат физико-математических наук, //Ц}0 З.М. Иванченко

/63

мтоо

Общая характеристика диссертации

Актуальность темы

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

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

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

Научная новизна

На основе нелинейных дифференциальных уравнений в частных производных предложена математическая модель для описания процессов тепловой релаксации с учетом плавления в условиях мощного энерговыделения, не опи-

сывающихся в рамках традиционной постановки задачи Стефана. Разработан новый комплекс алгоритмов и программ для исследования процессов треко-образования в высокотемпературных сверхпроводниках. Выполненные численные эксперименты впервые показали с большой степенью достоверности, что механизм температурного пика является ответственным за процессы трекооб-разования в иттриевых высокотемпературных сверхпроводниках. Разработана программа расчета термоупругой релаксации в образцах, подвергшихся облучению мощным пучком ионов, отвечающим источникам со взрывной ионной эмиссией. Численные эксперименты, выполненные в рамках используемой модели, показали, что механические напряжения, возникающие вследствие теплового расширения вещества, не являются ответственными за модификацию образца. Предложена модель высокотемпературной электрон-ионной релаксации в плазме, образующейся при схлопывании кавитационного пузырька в дейтерирован-ном ацетоне, и разработан програмный комплекс для нахождения неизвестых физических параметров модели. Расчет времени электрон-ионной релаксации подтвердил гипотезу существования перегрева ионной компоненты, принципиально важную с точки зрения возможности осуществления реакций термоядерного синтеза в процессах схлопывания кавитационного пузырька в дейтериро-ванном ацетоне.

Практическая ценность

Программы были использованы для описания экспериментальных данных по трекообразованию, полученных в нескольких ускорительных центрах (ОИЯИ, GANIL и др.), для объяснения процессов обработки поверхностей материалов интенсивными ионными пучками (ОИЯИ, ЛВЭ) и для объяснения выхода продуктов ядерных реакций в кавитационных пузырьках в дейтериро-ванном ацетоне (эксперимент в Окридже, США).

Апробация работы

Основные результаты диссертации докладывались на семинарах ЛИТ ОИЯИ, Дубна, были представлены и докладывались на международных конференциях: "International Scientific Conference on Mathematics", (г. Херляны, Словакия, 2000); "European Network on Ion Track Technology", CIRIL-GANIL, (г. Казн, Франция, 2002); "7th international Scientific Conference "(г. Кошице, Словакия, 2002); "The First International Congress on Mathematical Modelling", (г. Дубна, 2002); "V Международный уральский семинар по радиационной фи-

зике металлов и сплавов", (г. Снежинск, 2003); "XVI Всероссийская конференция по проблемам математики, информатики, физики и химии" (г. Москва, 2005).

Публикации и Личный вклад автора

По результатам диссертации опубликовано двенадцать работ. Результаты, выносимые на защиту, получены лично автором. В работах, опубликованных в соавторстве, личный вклад автора был определяющим на всех этапах выполнения данных работ.

Объем и структура диссертации

Диссертация состоит из введения, трех глав, заключения и приложения. Список литературы содержит 98 наименований. Полный объем диссертации -109 страниц машинописного текста, включая четыре таблицы и тридцать восем рисунков.

Содержание работы

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

В работе рассматривается ряд современных применений системы нелинейных дифференциальных уравнений в частных производных вида

ПЕ

Р-^ = -Ре? ■ v + V ■ (КеЧТе) ±9-(Те- Т.) + де(г,«), (1) Г) Ф

рт = 'у + у'{куТ{) т 9'{Те ~ т,) + (2)

где р, V, Р, Е, Т, К - плотность, скорость, давление, энергия, температура и теплопроводность, электронной (е) и ионной (г) подсистем, д - постоянная электрон-атомного взаимодействиея, д(г, £) внешний источник мощности,

= §1 + V ■ V - полная производная. Граничные и начальные условия формулируются в контексте конкретных задач.

В первой главе диссертации предпринята попытка построения реалистической модели термоупругих процессов в материалах, испытывающих интенсивную ионную бомбардировку. С этой целью развита модель пространственно-

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

В первом параграфе с целью изучения явлений формирования волн сжатия была сформулирована модель, описывающая совместную эволюцию термоупругих напряжений и температур в образце, на базе следующей системы уравнений с начальными и граничными условиями:

д2а . . д2Т , да(х,0) п

(1+А.Т)^ = (4)

„(О,*)-(М> = 0, ТМ) = 1, =

В уравнениях (3) - (4) все величины даны в безразмерном виде. Безразмерные величины можно выразить через физические постоянные параметры по формулам:

Л0т _ Ет2 ЕатТ0 _ ат<т0 _ Еа^Тр

А — 75, А1 — -¡х, А<± — , Аз — , Л4 — ,

СоРоЦ Ро'о °о СоРо соро

где До - теплопроводность металла, Со - удельная теплоемкость, ро - плотность, Е - модуль Юнга, ат - коэффициент объемного расширения. Функции Т = Т(х,Ь), а = <х(х,£) описывают температуру образца и напряжение, д = д{х,Ь) -мощность теплового источника.

В первом пункте первого параграфа сформулирована модель энерговыделения ионов с учетом конкретных особенностей импульсных источников ионов со взрывной эмиссией, а также точной пространственно - временной картины энергетических потерь пучка ионов в веществе.

Зависимость от времени интенсивности падающего пучка в некотором приближении (передающем, по меньшей мере, ее качественное поведение) можно аппроксимировать формулой

1 — е~'/т1

где Т\,Т2 ~ (1 - 5)-Ю-8 с - длительности нарастания и спадания импульса, At - его длительность. Аналогичным образом описывается и мощность, Р = 1^)1/, попадающая на поверхность образца: Р — Р0/(£), где II - ускоряющее напряжение. Здесь предполагается выполненным следующее очевидное условие самосогласованности: /•»

Р0 / /(О <и = <?0, ./о

для указанных интервалов значений плотности выделившейся мощности Р0 и плотности выделившейся энергии <5о-

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

= (5)

ах

В частности, для ионов углерода, сталкивающихся с железным образцом, в системе единиц [Е] = кэВ , [гс] — м, следует взять а — 5,82 • 107. Рассматривая соотношение (5) как уравнение, описывающее зависимость энергии налетающего иона с начальной энергией Е0 от пройденного пути, находим

ад = (¿Г - ¡*)2.

Время прихода иона в точку, отстоящую от поверхности на расстоянии х, равно

Гх ¿х \ Гх <1х 2

Ф)= -¡-г = т -г■=-- = --г1п(1-а;/ад>)), 0 < х < ВД),

Уо ЬУо ЕУ — ах/2 аЬ

где символом Я(Ео) обозначен пробег иона в веществе:

Г*> АЕ 2 ^

R(F\- f dE V

R{Eo)-J0

а коэффициент Ъ для ионов в Fe равен b = 1,22 • 106, если [г>] = м/с. Распределение по времени выделяющейся мощности внутри образца с учетом времени

движения иона до точки х можно теперь представить в виде о р

^(z, t) = щщ (1 - x/R(Eo)) f(t - т(®)) 9 (t - т(х)). (6)

Во втором пункте приводится метод численного решения. Система (3) - (4) решалась на прямоугольной сетке переменных х и t с шагами hx и т( в интервале х € [0,1], t € [0, tntax\. Длина х нормировалась на толщину железной мишени 1о, а время t - на длительность импульса г = 3 • Ю-7 с.

Т(х.ц

Т{х,1)

Рис. 1: Пространственно-временное распределение температуры в образце (в безразмерных единицах).

о(х.Ч

Рис. 2: Пространственно-временное распределение напряжения в интервале времени от 0 до Ю-8 с.

Для решения системы применялась следующая конечно-разностная схема с весами 71,72 € [0,1]: сг*+1 — 2(7* + сг*-1 + +

-^-= Лх[71-«-+ (1 - 271)-Л2-+

I -к-1

т1+1

-V 1"-Т2 ' (7)

т->к+1 т* плк+1 <урк+1 . грк+1 грк _ гурк ±_ тк

(1 + Л4Г*) ^ = ¿2 + + (1 - 72) %

¿к

-ЛзГ,^' (8)

п

где а,* = аОг.Л), 7* = Т(х,,Ьк) и <ь+1/2 = Я(х„1к + ъ/2).

В третьем пункте приведены результаты численных расчетов. На рис. 1 показано распределение температуры в металлическом образце толщиной 10~® см вплоть до момента времени 9 ■ Ю-7 с при энергии частиц падающего пучка 300 кэВ. На рис. 2 приведено распределение напряжения в образце при той же энергии. Отчетливо видна волна разрежения, сформировавшаяся в результате отражения волны сжатия от задней поверхности.

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

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

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

КаЫ^±М _ КьщШЕ|zM = LpsalVs. (9)

Здесь Vs = d£s/dt - скорость движения границы S, Ksoi и Kitq - коэффициенты теплопроводности материала для твердой и жидкой фаз, L и psoi - удельная теплота плавления (или энтальпия фазового перехода) и плотность соответственно. Полная математическая формулировка задачи Стефана включает, кроме (9), условие, учитывающее тот факт, что фазовый переход происходит при постоянной температуре,

T\s = Т, (10)

где Т* - температура плавления, а также - закон сохранения энергии:

рС = —div j + ç(x, t).

<?(x, t) описыват мощность внешнего источника тепла, а С - коэффициент теплоемкости.

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

(рС + L 6(Т — T*)) + v grad Tj = div(AT grad T) + 9(x, t),

где Ь 6{Т — Т*) дТ/дЬ дополнительный вклад в теплоемкость тепла, израсходованного на фазовый переход, V grad Т учитывает возможное температурное изменение, обусловленное конвекцией (в задачах, рассматриваемых в данной главе, не возникает необходимости учитывать этот член). Таким образом, главная

идея этого подхода сводится к предложению учета удельной теплоты плавления Ь в качестве дополнительной компоненты теплоемкости рС, которая, однако, дает вклад только в точке фазового перехода, когда Т = Т*. Но даже в этом подходе уравнение (11) обычно считается только следствием условия (9).

В первом пункте показано, что на самом деле условие (11) дает правильное описание фазовых переходов даже в том случае, когда соотношения (9) и (10) не применимы.

Во втором пункте рассмотренны фазовые переходы в присутствии мощных источников энергии. В частности, изучалось уравнение теплопроводности

= тЛк{Т)^)+ч (12)

с начальными и граничными условиями Т(х,0) = То, = — 0. Чис-

ленный эксперимент проводился на мишени из железа. Временная и пространственная зависимость мощности источника (в безразмерных единицах) была взята в виде:

д(х,Ь) = (3 д1{х)д2{г), где дг(г) = 1

1 + ехрр,(г - г,)'

Здесь <3 описывает полную энергию источника, ф = 59,44, Х\ = 0,07, = 1, /1, = 100. Небольшое различие между физическими параметрами для твердой и жидкой фазы в (12) не учитывалось в виду его несущественности для целей данного исследования.

Уравнение (12) решалось численно на прямоугольной пространственно-временной сетке переменных х и Ь с постоянными шагами /гт и соответственно. Для решения применялась конечно-разностная схема с весами для уравнений параболического типа:

грк+\ _ Фк

1 « 7; и

>4

(13)

где I и

7 е [о, 1], т} = т(х„1к), е* = р(Т*)с(Т*), д*+1 = д(х},1к + -£),

На рис. 3 отчетливо видно формирование "плато", высота которого соответствует температуре плавления, для пространственного распределения температуры. Узкая полоса, ограниченная двумя разрывными линиями на рис. 3 и 4, отвечает ширине сглаживания ¿-функции, введенного при численном решении уравнения (И). Наличие двух скачков пространственной производной температуры несколько замаскировано именно этой размазкой ¿-функции. Рис. 3 демонстрирует временную зависимость температуры на двух разных глубинах

Рис. 3: Зависимость температуры от х в моменты времени: 1) £ — 0,16; 2) Ь = 0,18; 3) < = 0,20; 4) I = 0,21 (Д = 0,01).

Рис. 4: Зависимость температуры от 4 на глубинах: 1) х 0; 2) х = 0,04 (Д = 0,01).

образца. Очевидно, что такое поведение температуры существенно отличается от традиционного описания фазового перехода в рамках классической задачи Стефана (9) и (10).

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

В первом параграфе приведена постановка задачи. Тепловые потоки в электронной и атомной подсистемах в рамках модели температурного пика описываются с помощью следующей системы связанных нелинейных дифференциальных уравнений:

Же

ЯГе \д рСе(1е)— = - — т гаг

гКе(Те)

дг

^■>141

дТ

-0-(Те-ТО + д(г,О, (14)

+ д.(Те-Т,), (15)

где г радиус в цилиндрических координатах, Те и Тг - температура электронов и решетки, Сс,, и Ке>г - теплоемкость и теплопроводность электронной и атомной подсистемы, д - постоянная электронно-атомного взаимодействиея, д(г, ¿) - источник мощности электронной подсистемы.

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

Уравнения (14), (15) не учитывают г зависимость Те и Т„ так как изменения потерь энергии иона по г малы. Это видно, например, из того, что величина пробега иона в веществе до его полной остановки существенно превышает радиус трека. Это обстоятельство позволяет ограничиться рассмотрением задачи в двумерном пространстве, а учет радиальной симметрии, имеющей место тогда, когда облучение ведется перпендикулярно плоскости (001) кристалла УВагСизОу-х, сводит задачу к одномерному случаю в цилиндрической системе координат.

Основные уравнения модели (14) - (15) фактически представляют собой закон сохранения энергии, записанный по отдельности для электронов и атомов. В (14) энергия, получаемая электронами при прохождении иона через вещество, учтена с помощью внешнего источника д(г, I). Атомы, в свою очередь, согласно (15), нагреваются от более горячих электронов за счет электронно-атомного взаимодействия, представленного слагаемым у ■ (Те — Т,). Соответсвующие потери энергии электронами учтены в (14) с помощью такого же члена, взятого с обратным знаком.

Начальные и граничные условия были взяты в виде:

где То температура окружающей среды, параметр обрезания rmm = 0,1 нм вводим из-за сложности вычисления источника электронной подсистемы в точке г = 0 (да и само понятие температуры на масштабах, меньших размеров атома, в рамках рассматриваемой модели не вполне корректно). Параметр rmax = ЮО, 1 нм взят в качестве физической бесконечности. Численные экс-

Гв(г,0) = Г,(г,0) = ГО)

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

Во втором параграфе описана дельта-электронная модель энерговыделения. Ее дальнейшее развитие было выполнено в направлении учета временной динамики процесса диссипации энергии, запасенной в ¿-электронах. Время, за которое электрон достигает точку на расстоянии Ъ от центра траектории иона, равно

т ~ /о v(b) - Л_ь у(П -г)~ с Л(д_ь) ЛЕ {м) [Е(Е + 2тес2)]1/2' (16) где Д = Ь + г (г - остаточный пробег ¿-электрона с энергией Е в точке Ъ), г = г(Е) соотношение пробег-энергия для электронов в материале. Плотность энергии выделившейся в момент времени I в объеме (2л- Ъ йЪ х единица длины) определяется формулой

где Е{Ь,Ь) - решение уравнения (16), йИ/йЕ - число ¿-электронов на единицу энергии, рассчитанное с помощью формулы Резерфорда. Зависимость пробега от энергии г(Е), а также обратная зависимость Е(г) были получены путем аппроксимации известных экспериментальных и теоретических данных.

На рис. 5 приведена плотность энерговыделения в электронной подсистеме УВагСизС^-я ионом 129Хе(2,6 МэВ/нуклон). Около 80 % полной энергии <1Е/<1х выделяется к моменту времени Ю-15 с, причем большая ее часть в области г < 1 нм при < < 0,15 фс, хотя процесс энерговыделения продолжается вплоть до I ~ Ю-5 с и г ~ Ю-3 см.

Тормозная способность, рассчитанная в рамках данной модели как радиальный интеграл от распределения дозы, находится в соответствии со значениями, найденными по программе БММ 2003 в пределах 12% точности. Хотя такая несогласованность отражает реальные возможности теории, существующей в этой области, радиальное

Рис. 5: Плотность энерговыделения в электронной подсистеме.

распределение энерговыделения было перенормировано к значениям БММ, которые обычно принимаются в качестве эталона.

В третьем параграфе приводится описание электронной системы. Как уже упоминалось выше, основные уравнения модели (14) - (15) являются законами сохранения энергии. Они допускают, как квантовую, так и классическую интерпретацию, фактически содержащуюся в теплофизических константах, в частности - теплоемкости и теплопроводности электронной и атомной подсистем. Теплоемкость электронов в широком температурном интервале может быть найдена численным методом по формуле

рСе(ТЕ) = Jе dn(e),

где /(г, Те) - распределение Ферми, dn(e) = т](е) de, and 77(e) - плотность уровней электронов в УВа2Сиз07_х. Таким образом полученное значение параметра Зоммерфельда, 7 = рСе/Те, с учетом имеющихся экспериментальных данных для плотности электронных уровней, равно 2,4±0,8-10-4 Дж/(см3-К2).

Зависимость параметров рассматриваемой модели от температуры электронов можно учесть, если принять во внимание теорию Аллена и экспериментальные данные по времени электронно-атомной релаксации в YBa2Cu307_:r. В самом деле, теплопроводность электронов связана с параметром температуропроводности по формуле

Кс = DepCc, где рСе=-уТе, 7 «2,4- Ю-4 Дж/см3К2.

Эффективное время электронно-атомной релаксации г равно

т = pCJg.

Поскольку Се является линейной функцией от Те, функция т(Те) имеет в рассматриваемом случае такую же линейную форму г = (7/ff) Те = аТе, как это было предсказано в теории Аллена, где

_ тг kp т Т 3 А' < w2 >

Используя экспериментальное значение А' < и2 >= 475 ± 30 меВ2, можно определить параметр а:

а=( 1,28 ±0,06)-10"16 с/К.

Таким образом, неизвестную постоянную электрон-атомного взаимодействия д в уравнении (14) можно выразить через параметры а и 7:

д = ч/а-

С вычислительной точки зрения наиболее удобным оказалось представление уравнения (14) для электронной компонеты системы в виде:

рС^Ж - гд~г

rDcPCc{Tc

&j&.{Te-Tt)+q(r,t). (17)

В четверток параграфе описанна атомная подсистема и ее параметры. Температура плавления Тт в УВагСизОт-г, найденная in situ методами нейтронно-дифракционного анализа, приблизительно равна 1070 °С.

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

YBa2Cu307-i —» solid + liquid + О2

от давления кислорода. Обычное рассмотрение, основанное на уравнении Кла-пейрона-Клаузиуса, позволяет найти изменение энтальпии, отвечающее данной реакции: Qm = 810 ± 5 кДж/мол.

Для теплоемкости решетки, согласно закону Дюлонга и Пти, было взято значение рСг = 3,1 Дж см-3 К-1, где р = 6,39 г см-3. Теплопроводность атомной системы Кг была выбрана в соответствии с существующими работами приве-деными в дисертации, К, — 5,6 • Ю-2 Дж(с см К)-1. Поскольку величина К, считается не зависящей от температуры, можно ввести коэффициент температуропроводности Dx — К,/рС, и уравнение (15) записать в виде:

f = Д>"(Г,)ДГ. + _ Tt)> (18)

где

С,е// = С, + Qm 8(Тт - Ti) (19)

— эффективная теплоемкость, которая включает теплоту плавления Qm — 1,216 кДж/г. Формально в (18) можно положить (с учетом того, что 5 функция находится в знаменателе выражения для D,)

о:"(т{ ф Тт) = А, ю:п(тг = Тт) = о

и, таким образом, представить (18) в виде "прямой суммы" двух простых уравнений теплопроводности с регулярными коэффициентами.

Для численного решения системы (17) - (18) было проведено сглаживание коэффициентов С'^ и О'^ = К,/рСв некоторой окрестности температуры плавления (в интервале Тт — Д < Т, < Тт + Д) следующим образом:

(20)

с параметрами а = 5 К, Д = 4,5 <т.

В пятом и шестом параграфах приведен метод численного решения и проверки точности разностных схем. Система (17) - (18) решалась на неравномерной прямоугольной сетке с постоянными шагами h\, h^ по радиусу г € [г™,, Гтах] И ПОСТОЯННЫМ ШвГОМ ПО времени i € [0, ¿щах].*

¿А = {г, = rmm + г • hui = Q,...,Ni,r, = rmld + i ■ h2,i = iVj + 1,..., N2} йы = {<• = i-ht,i = 0,...,N}

hl = {rmid - rmin)/N\, h2 — (Гтах - rmtd)/N2, ht = tniax/N По нескольким причинам мы поставим rmui = 10,1 нм и возьмем hi < h2. Это, во-первых связано с тем, что большая часть энергии источника выделяется в области г < 10 нм. Во-вторых, граница фазового перехода, которую мы рассматриваем находится также в этом диапазоне. Таким образом удается более точно и экономно решать систему (14) - (15). Введем обозначение:

Се = рСе(Тс), Ке = DepCe(Te), = pc:ff, А = Ас//.

Для аппроксимации системы уравнений использовалась конечно-разностная схема с весами 7 £ [0,1]:

т_ Tfc 1 С k

С,4 ht 6J = ¿А + (1

(21)

rpk+l _ npfc ri /-< fc

lj , '' = [7 Tj T^ + (1 - 7) r, T,*] + -¿^-(7^ - T,*), (22) nt ~3 T? G,J

где

А(г,ВД) =

h3

А, = АЬ J = 0,...,A^! - 1, hj =h2,j = Nu...,N2, fr} — (h} = h3+i)/2,

T*=T(rj,th), Ck3=C(T*), K) = K(Tj), т* = t(Tj), q?1'2 = q{rj<tk+ht/2),

+ - K(Tj) + K(Tj+1) rj+i =-ö-' -2-'

Рис. 6: Т%{г,Ь) для иона шХе (2,6 МэВ/нуклон) в УВагСизОг-я. Разрывная линия отвечает Т, = Тт, непрерывная линия соответствует Т, = Тт + Д, т.е. полному окончанию процесса плавления.

Рис. 7: Бифуркация траекторий описывающих электронную температуру в точках, где Т, = Тт + Д, Т4 = Тт и Г, = Тт — Д (кривые о, Ь и с соответсвенно) в УВа2Си307-1 для 129Хе (41 МэВ/нуклон).

В седьмом параграфе описаны особенности найденного решения. Численное исследование системы (17) - (18) с физическими параметрами, соответствующими УВа2Сиз07_х, выявило ряд новых, достаточно неожиданных эффектов. Один из них, о котором уже шла речь в предыдущей главе диссертации, обусловлен присутствием источника мощности в исходных уравнениях.

На рис. 6 показано распределение температуры Т,(г, £) для иона 129Хе с энергией 2,6 МэВ/нуклон. Вместо отчетливой границы двух фаз, видно пространственный слой ненулевой толщины, который нагрет до температуры плавления, и два скачка градиента температуры на внутренних и внешних границах слоя.

Численные эксперименты также выявили пороговое явление, имеющее место, когда значения Юе п Qf являются достаточно большими (при фиксированной величине параметра £>,). Этот случай показан на рис. 7, где приведены температуры электронов в точках г(<), отвечающих следующим фиксированным температурам атомов: Т;(г,¿) = Тт + А, Тг (г, ¿) = Тт и Т,(г,<) = Тт - А (кривые а, Ь и с, соответственно). Точка на кривой а отвечает теоретическому радиусу трека. Видно, что при небольших изменениях Ое электронные температуры вдоль траекторий с и 6 испытывают бифуркацию. В результате, это малое изменение пареметра Д. приводит к резкому изменению траектории а, описывающей формирование трека.

Причины такого нерегулярного поведения разъясняют рис. 8 и 9, где показа-

Рис. 8: Те(г), Г;(г) для иона 208РЬ с энер- Рис. 9: Те(г), Т,(г) для иона 129Хе с энергией 3,7 МэВ/нуклон в УВа2Сиз07_х. гией 10 МэВ/нуклон в УВа2Сиз07_х.

но распределение электронной температуры и температуры решетки в момент времени Ьа, когда радиус расплавленной области достигает своего максимального значения г = а (обозначается граница трека). Радиус трека обозначен точкой на непрерывной кривой. "Плато" при Т, = Тт свидетельствует о существовании в этот момент жидкой фазы вещества, возникшей в результате электронного нагрева. Однако, к моменту времени Ьа электронная температура в центре трека уже стала ниже чем Тш, и, таким образом, начался процесс "электронной закалки" материала (см. рис. 8). В некотором смысле противоположный процесс, тогда атомы в момент времени ¿а все еще нагреваются электронами, который возможен при других экспериментальных условиях, показан на рис. 9. Численные эксперименты показали, что переход от электронного нагревания к электронному охлаждению, вызванному малыми изменениями £>е, и является причиной неустойчивости, показаной на рис. 7.

В восьмом параграфе проведено сравнение расчетов с экспериментом. Экспериментально полученные радиусы треков гехр, в УВагСизОг-* с ориентацией оси [001] параллельно траектории иона приведены в таблице 1 вместе с результатами наших расчетов. Температуропроводность электронов Ц. = Ке/рСе выбрана таким образом, чтобы теоретические значения радиусов треков совпадали с экспериментальными гехр. Неточность расчетов для 120Хе (41 МэВ/нуклон) обусловлена его попаданием в область бифуркации (см. рис 7). *

Расчетные зависимости температуры решетки в области трека показаны в виде двух нижних кривых на рис. 10, где пунктирная и непрерывная линии относятся к границе (г = радиусу трека) и центру трека (г = 0) соответственно.

Таблица 1: Экспериментально найденные радиусы треков ГехР в монокристале УВа2СизС>7_х и результаты теоретического описания. Выделившаяся энергия ¿Е/скс найдена с помощью программы ЭШМ 3003.

Ион Энергия, аЕ/ах, Техр-> а, Де,

МэВ/нуклон кэВ/нм нм нм см2/с

129Хе 1,3 26,2 2-3 2,71 0,730

129Хе 2,6 30 2,5 2,49 0,768

129Хе 10 27,9 1,3 1,35 0,605

129Хе 27 18,7 1,3 1,6 0,326

129Хе 41 14.8 0,56 0,44-1.55 0.223-0,222

208рЬ 3,7 43,7 4 4,1 1.130

208рЬ 10 42,5 3 3,02 1,015

208рЬ 20 37 3,5 3,52 0,805

208рЬ 25 34,5 3 3,06 0,732

* [Ю-11 с]

Рис. 10: Эволюция температуры решетки на границе и в центре река (пунктирная и непрерывная линии соответственно) для иона 129Хе (2,6 МэВ/нуклон).

04-

«Л

/

01-1-i-1-1-i-

0 1 2 3 4 Ъ

Те [104 К]

Рис. 11: (а) Зависимость электронной температуропроводности Бе от <1Е/йх в УВагСизОт-ц найденная с помощю МТП (точки). (Ь) Теоретическая зависимость Д.(Те) для аморфного углерода.

Влияние электронной закалки на скорость охлаждения очевидно из сравнения этих кривых с двумя верхними, полученными для случая, когда обратная передачи тепла от атомов решетки к электронам, имеющая место при Те <Тг, была искусственно "выключена".

Рис. 10 позволяет проверить гипотезу "эпитаксиального восстановления", согласно которой внешняя часть трека не переходит в аморфное состояние из-за краткости времени тт ее пребывания в жидком состоянии. Значения для времен пребывания при температуре выше Тт + Д действительно сильно различаются (точки а, на рис. 10). Однако, соответствующие значения тт для точек, в которых Тт и Tt = Тт — А очень близки (точки Ь, и с„ соответственно). Из выражения (20) видно, что теплота плавления поглощается главным образом в узком температурном интервале Тт — а < Тг < Тт + а в окрестности температуры плавления Тг = Тт, поэтому разумно считать, что в рассматриваемой модели для всех точек трека rm ~ 0,4 • Ю-11 с (см. продолжительность времени пребывания в расплавленном состоянии для точки Ь3, отвечающей Тг = Тт). Таким образом вычисления показали, что очень существенное "эпитаксиаль-ное восстановление", обусловленное малым временем тт пребывания внешних областей трека в расплавленном состоянии, в рамках традиционной модели термопика представляется маловероятной.

Зависимость полученной электронной температуропроводности De от потерь энергии dE/dx налетающего иона показана на рис. 11 (а). Непрерывная линия представляет соответствующие усредненные значения. Разрывная линия демонстрирует неточность, обусловленную экспериментальными ошибками параметра а = (1,3 ± 0,1) • 10"16 с/К. Экспериментальная ошибка теплоты плавления незначительна (представлена маленькими "волнами" вдоль непрерывной линии). Основными источночниками ошибок, наблюдаемых в виде точек, рассеянных вокруг сплошной линии, являются флуктуации экспериментально измеренных радиусов треков, а также ошибки расчета величины dE/dx.

Рис. 11 (а) свидетельствует о том, что параметр Д, для YBa¿Cu307_x нельзя считать независимым от электронной температуры, как это предполагается в каэнской версии МТП. В поддержку этого вывода на рис. 11 (Ь) приведена зависимость Д. от Тс для аморфного углерода, полученная на основе теоретических результатов работы [G. Sciwietz et ai, Nucí. Instr. and Meth. B164-165, 354 (2000)]. Видно, что в этом случае температуропроводность также увеличивается с ростом температуры электронов в области Те ~ 103 К. Именно такие температуры являются типичными при трекообразовании в YBaíCu307_:t.

Наиболее важно, однако, то, что найденные значения Д, оказались действительно близкими к величине Д. ~ 1 см2/с, которая обычно считается вполне обоснованной теоретически в каэнской и некоторых других моделях формирования треков. Таким образом, результаты расчетов, полученные в рамках разработанной в диссертации модели температурного пика, выглядят с этой точки зрения достаточно реалистичными.

В третьей главе проведена проверка гипотезы перегрева, а также связанное с этим развитие модели в направлении явного учета процессов охлаждения ядер электронами в сверхплотной сильно неравновесной плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне. Актуальность этих исследований обусловлена появившимися в последнее время сообщениями о регистрации продуктов термоядерных реакций при акустической кавитации в C3DeO [Taleyarkhan R.P. et al., Phys. Rev. E. 2004. V.69. P. 036109]. Поскольку в рассматриваемой модели предполагается, что электроны характеризуются определенной температурой Те, расчеты энергии электронов Ее, фактически сводятся к расчетам их теплоемкости Се{Те). Информация о параметре электрон-ионной связи д, как и в процессе трекообразования, эквивалента знанию времени т электрон-ионной релаксации. Основные проблемы, связанные с предлагаемым в диссертации обобщением модели [Nigmatulin R.I. et al., J. Acoust. Soc. Am. 2003. V.113. P. 2205] связаны именно с необходимостью расчета параметров g и с(Те), так как информация об остальных параметрах уже содержится в той или иной форме в работал других авторов.

Первый параграф посвящен оценке времени охлаждения ядер электронами в сверхплотной сильно неравновесной плазме.

В первом и втором пунктах приведена оценка начальной температуры и дано краткое описание общей теории процессов электрон-ионной релаксации в высокотемпературной плазме.

В третьем и четвертом пунктах дана конечно-разностная аппроксимация, традиционно используемая при расчетах электрон-ионной релаксации, с учетом зависимости от Те. Для простой оценки времени релаксации может быть использовано уравнение Ландау:

Ç = Те*~Те. (23)

at Tjve

Предполагая, что температура ядерной подсистемы в начале процесса релаксации равна 2 • 108 К, получим для плазмы, образовавшейся в результате ионизации молекул C3D60, сотношение между температурами электронной и ядерной

компонент:

Ты =2 - 108-3,2Те.

Условие равновесия Тщ — Те = Тсд дает Т^ = 10®/2,1 К. Учитывая тот факт, что время релаксации внутри ядерной подсистемы существенно меньше времени ядерно-электронной диссипации, получаем для характерного времени

охлаждения ядер электронами соотношение

1 _1_ + +

Ты с То тс То

Оценка величины г получена в [Спитцер Л., Физика полностью ионизированного газа. М.: ИЛ, 1957.1: ,,,

250АДГ

~~ пг2% 1п Л,

Здесь А, - атомный вес ядра, 2г - его атомный номер, 1пА, - кулоновский логарифм, который может быть вычислен лишь приближенно (разные авторы приводят разные значения). Если в качестве радиуса экранирования электрического поля в плазме взять радиус Дебая, то придем к следующей часто используемой оценке (в СГСЭ):

3 (кТе)3/2

2 ¿¿{Ьту-^Я!'

(25)

—• т.,

00 о

Е?

Для решения уравнения Ландау применялись схемы Рунге-Кутта, четвертого порядка точности. На рис. 12 показаны зависимости 7лг(£) и Те{Ь) К, отвечающие уравнению (23) при наг чальной температуре электронов Ге(0) = 4-10® К. На оси Ь отмечен момент времени, когда наг чальная температура ядер уменьшается приблизительно в е сы 2,72 раз. Видно, что неопределенность в начальной температуре электронов Те = 4 • 10® Ч- 107 К почти не сказывается на оценке длительности ядерно-электронной релаксации, которая оказывается приблизительно равной 8 • 10"13 с.

В пятом пункте сформулированны основные выводы. Оценки с использованием аналитической формулы (24), полученной Спитцером, показывают, что

1 2

1-12 ,

t [Ю-" с]

Рис. 12: Временная эволюция температуры ядер Т/у(£) и электронов Те(г) с учетом зависимости времени релаксации туус от Те(4)

электронная компонента плазмы действительно может оставаться относительно холодной на протяжении заметной части времени существования сверхплотного состояния вещества в схлопывакяцемся кавитационном пузырьке в СзБбО. Расчеты, выполненные на основе такого модельного допущения, предсказывают в этом случае возможность осуществления термоядерных реакций D-D слияния со скоростью порядка одного события на одно схлопывание пузырька при температуре D-ацетона около 273 К. Оптимизма относительного этого новой возможности реализации термоядерного синтеза прибавляют также два дополнительных обстоятельства. Во-первых, более тщательное численное моделирование показывает, что оценки, выполненные с помощью аналитических формул типа (24), использовашихся в данной работе, дают некоторое занижение величины времени релаксации. Во-вторых, расчеты в рамках гидродинамической модели, говорят о том, что дальнейшее понижение температуры D-ацетона приводит к уменьшению времени существования термоядерной плазмы почти до значения Ю-13 с, в то время как выход нейтронов при этом увеличивается.

Во втором параграфе выполнены расчеты теплоемкости электронов в неравновесной плазме, которые использовались при оценке температуры электронов в момент начала электрон-ионной релаксации.

В первом пункте приведена формулировка модели для вычисления теплоемкости электронов в плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне и приведены результаты расчета.

Среднее число электронов с энергией Е будем описывать распределением Ферми для невзаимодействующих частиц:

где Т - температура электронов, К в - постоянная Больцмана, ц(Т) - химический потенциал электронов, находящихся в электрическом поле ионов.

Из того, что при Т = 0 все электроны находятся в связанном состоянии следует, что значение химического потенциала ц при нулевой температуре равно наибольшей энергии связи электрона в атомах С, Б и О, т.е. 871,1 эВ. Выделяя в выражении для химического потенциала эту величину отдельным слагаемым,

(N) =

1

(26)

е(Е-„(Т))/КвТ + I

Г

м(Т) = М0) + Д/|(П

выражению (1) нетрудно придать вид

s

1

e(tk-bn)/KBT + I

соответственно для связанных и свободных состояний. Здесь е, = /л(0) — Е, ~ энергия связи, отвечающая г—му электронному уровню, е* = р\/2т - энергия свободного электрона.

Зависящая от температуры поправка к постоянной части химического потенциала Др.(Т) находится из условия сохранения полного числа электронов

- <ЛГ,,ь)) = (27)

г к

где суммирование выполняется с учетом кратности вхождения атома данного сорта в молекулу С3Б60. Заменяя обычным образом сумму по состояниям свободных электронов интегралом (здесь V - объем, приходящийся на одну молекулу вещества)

Уйут3^ \ftdt _ у/2у{тКвТ)3'2 ф<1г

тг2/г3 Уо е('"дМКвТ + 1 " п2цз у0 ег-^/квт + ! '

получим уравнение для определения величины поправки Ар.(Т)

Е^О- + ])= СОП8{ГЗ/2/°° + ! . (И)

где п, - кратность вхождения атомов в формулу вещества.

Теплоемкость электронов, отнесенная к одной молекуле вещества, может быть найдена прямым дифференцированием полной энергии свободных электронов: . . ^ .

где суммирование выполняется по всем связанным состояния электронов.

На рис. 13 представлена рассчитанная таким образом теплоемкость электронов в области сильной ионизации при сжатиях р/ро =1, 10 и 100, переведенная в более удобные для практического использования единицы. Видно, что она не описывается линейной зависимостью С = 7Т, справедливой для вырожденного электронного газа в широкой потенциальной яме с плоским дном при температурах меньше или порядка температуры Ферми, когда постоянная (Зоммер-фельда) 7 может быть рассчитана по формуле

г-2 /ЯЛ2/3 тп „ ( V \2/3

где N - полное число электронов в молекуле.

Для моделирования процессов схлопывания кавитационого пузырька важно иметь простые, но в тоже время, достаточно точные интерполяционные формулы расчета зависимости Су(Те,р/ро).

Рис. 13: Теплоемкость электронов для молекулы D-ацетона при разных сжатиях р/ро, рассчитанная в соответствии с рассматриваемой моделью.

Во втором пункте решена задача интерполяции с точностью, достаточной для математического моделирования процесса. Найдены интерполяционные формулы, которые позволяют с помощю нескольких вычисленных значений теплоемкостей при сжатиях р/ро =1, 10 и 100 находить значения во всем интервале температур и плотностей образующейся плазмы, с относительной точ-ностю меньше 4%.

В заключении сформулированы основные результаты диссертации.

В приложении к диссертации находятся следующие материалы:

• Таблица значений химических потенциалов для некоторых характерных температур электронов при разных степенях сжатия, которые позволяют расчитать, используя распределение Ферми, средние числа заполнения электронных уровней в D-ацетоне.

• Программа 1, которая вычисляет значения теплоемкостей электронов в СзОеО решая уравнение (28). Она легко может быть адаптирована для проведения аналогичных расчетов и для других веществ.

• Программа 2, которая вычисляет значения теплоемкостей электронов > в C3D60 с помощью интерполяционных формул.

Программы написаны в среде SCILAB 3.0-RC1 (May 17, 2004), которая является 1 бесплатно распространяем аналогом MATLAB (http://scilabsoft.inria.fr).

На защиту выдвигаются следующие результаты

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

2. Разработан комплекс математических моделей и программ для расчета процессов тепловой релаксации в образцах, подвергаемых облучению ускоренными ионами. С его помощью были выполнены исследования взаимодействий ионов с веществом и получены следующие результаты:

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

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

в) В рамках концепции температурного пика рассмотрены процессы тре-кообразования в высокотемпературных сверхпроводниках. С учетом имеющихся теоретических и экспериментальных данных разработана модель процесса энеговыделения, описывающая его пространственно-временную динамику. Уточнена модель электронной подсистемы с использованием теории Аллена и результатов фемтосекундных лазерных экспериментов.

г) Продемонстрирована применимость концепции температурного пика для описания процессов трекообразования в иттриевых сверхпроводниках. Установлена бифуркационная зависимость решений от параметра температуропроводности электронов, связанная с процессом электронной закалки.

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

По теме диссертации опубликованы следующие работы:

1. Ayrjan Е.А., Fedorov A.V., Kostenko B.F., PribiS J. Calculations of temperature fields in the vicinity of ion track within a thermal-spike model. Journal of Computational Methods in Sciences and Engineering (JCMSE), 2,1-2, 2002, pp. 163-168

2. Kostenko B.F., PribiS J., Puzynin I.V. Stefan Problem and Beyond, e-print: math-ph/0302044, 2003, to be published in Journal of Computational Methods in Sciences and Engineering (JCMSE).

3. Костенко Б.Ф., Прибиш Я. Математическое моделирование трекообра-зования в высокотемпературных сверхпроводниках. Вестник РУДН, сер. Прикладная математика, 2005, т. 4, No 1, с.75-87.

4. Kostenko B.F., PribiS J., Goncharov I.N. Thermal spike model of track formation in YBo^CvqOt-x. Preprint JINR E17-2005-61, Particle and Nuclei Letters, 1(130), Vol. 3, 2006, pp. 31-44.

5. Амирханов И.В., Айрян Э.А., Федоров А.В., Костенко Б.Ф., Прибиш Я., Сархадов И. Численное моделирование термоакустических процессов, генерируемых интенсивными ионными пучками в тонких образцах. Сообщение ОИЯИ, Р11-2000-271, 2000.

6. Костенко Б.Ф., Прибиш Я. Оценка времени охлаждения ядер электронами в сверхплотной сильно неравновесной плазме. Сообщение ОИЯИ, Р4-2004-42, 2004.

7. Костенко Б.Ф., Прибиш Я. Теплоемкость электронов в плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне. Сообщение ОИЯИ, Р11-2004-193, 2004.

8. PribiS, J. Numerical solution of system of partial differential equations which describes thermoelastic effects during irradiate. Proceedings of International Scientific Conference on Mathematics, 2000, Herl'any, Slovakia, pp.156-158, (Univ. Technol. KoSice, 2000).

9. B.F. Kostenko, J. Pribis, I.V. Puzynin, V.A. Skuratov, S. Zinkle, Numerical Solution of Heat Relaxation Processes within the Thermal Spike Model. Proceeding of Workshop "European Network on Ion Track Technology", CIRIL -GANIL, Caen, FYance, 24-26 February, 2002, p. 39.

10. Ayrjan E. A., Kostenko B.F., PribiS J. Numerical simulation of heat relaxation processes within thermal spike model. Proceedings of the 7th International Scientific Conference, Section Applied Mathematics, 2002, KoSice, Slovakia, pp. 16-20.

11. И.Н.Гончаров, Б.Ф. Костенко, Я.Прибиш., Расчет поперечных размеров треков в монокристалле YBoq Си3 От-х с учетом зависимости скорости передачи энергии решетке от температуры возбужденных электоронов. V Международный уральский семинара по радиационной физике металлов и сплавов, Снежинск, 2003, Аннотации докладов, с. 19.

12. Kostenko B.F., Pribis J., Математическое моделирование гпрекообразова-ния в высокотемпературных сверхпроводниках. XLI Всероссийская конференция по проблемам математики, информатики, физики и химии, 1822 апреля 2005, Тезис докладов - Физические секции, Москва, РУДН, 2005, с.36.

Получено 22 сентября 2005 г

i I

f f

I

• Шбб

РНБ Русский фонд

2006-4 16388

Отпечатано методом прямого репродуцирования с оригинала, предоставленного автором.

Макет Я. А. Киселевой

Подписано в печать 22.09.2005. Формат 60 х 90/16. Бумага офсетная. Печать офсетная. Усл. печ. л. 1,62. Уч.-изд. л. 1,75. Тираж 100 экз. Заказ № 55020.

Издательский отдел Объединенного института ядерных исследований 141980, г. Дубна, Московская обл., ул. Жолио-Кюри, 6. E-mail: publish@pds.jinr.ru www.jinr.ru/publish/

Оглавление автор диссертации — кандидата физико-математических наук Прибиш, Ян

Введение.

Глава 1 Численное моделирование термоакустических процессов, генерируемых интенсивными ионными пучками в тонких образцах

1.1 Моделирование процессов формирования и релаксации термоупругих напряжений

1.1.1 Модель энерговыделения

1.1.2 Метод численного решения.

1.1.3 Результаты.

1.2 Моделирование фазовых переходов в тонких образцах облучаемых ионными пучками.

1.2.1 Невозможность постановки классической задачи Стефана

1.2.2 Численное моделирование фазовых переходов в присутствии мощных источников энергии

1.2.3 Математическое моделирование трекообразования

Глава 2 Количественное описание трекообразования в УВагСизОу-^

2.1 Постановка задачи.

2.2 Модель энерговыделения.

2.3 Описание электронной подсистемы

2.4 Описание атомной подсистемы.

2.5 Метод численного решения.

2.6 Проверка точности разностной схемы.

2.7 Особенности решения

2.8 Описание радиусов треков

2.9 Выводы.

Глава 3 Высокотемпературная релаксация в плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне

3.1 Оценка времени охлаждения ядер электронами в сверхплотной сильно неравновесной плазме.

3.1.1 Оценка начальной температуры электронов.

3.1.2 Общая теория.

3.1.3 Конечно-разностная аппроксимация.

3.1.4 Время охлаждения ядерной подсистемы

3.1.5 Выводы

3.2 Теплоемкость электронов в плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне

3.2.1 Фомулировка модели и результаты расчета

3.2.2 Интерполяционные формулы.

Введение 2005 год, диссертация по информатике, вычислительной технике и управлению, Прибиш, Ян

Диссертация посвящена математическому моделированию динамики процессов высокотемпературной релаксации в веществе при облучении мощными ионными пучками. В работе рассматривается ряд актуальных применений системы нелинейных дифференциальных уравнений в частных производных вида

DE р-^ = -PeV • V + V • (KeVTe) ± д, {Те - ТО + qe(r, 0, (v.l) DE -P<V • v + V • (KiVTi) T9-{Te~ Ti) 4- ®(г, £), (v.2) где p, v, P, E, T, К - плотность, скорость, давление, энергия, температура и теплопроводность, электронной (е) и ионной (г) подсистем, д - постоянная электрон-атомного взаимодействиея, q(r,t) внешний источник мощности, ^ = §i + v - V - полная производная. Граничные и начальные условия формулируются в контексте конкретных задач.

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

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

Требования к ускорителям, предназначенным для генерации необходимых ионных пучков, определяются как поставленными технологическими задачами, так и особенностями взаимодействия ионов с материалами. В частности, для исключения активации мишеней необходимо ограничить энергию ионов величиной, находящейся ниже порога ядерных реакций. Для протонов, которые обычно всегда присутствуют в ионном пучке из-за различных загрязнений, это соответствует Е < 400 кэВ. Далее, важно также обеспечить существенную плотность выделяющейся энергии. Например, эксперименты показывают, что для того чтобы воздействие пучка на механические и антикоррозионные свойства поверхности металла стало заметным, необходимо, чтобы энерговклад пучка превышал энергию сублимации [2, 3], а для уменьшения влияния процессов тепловой релаксации важно добиться еще и достаточно высокой мощности поступающей энергии. Оценки показывают, что для практических целей обработки металлов часто бывает достаточно обеспечить длительность ионного импульса At ~ 10~7 с при энергии ионов Е ~ 300 кэВ и потоке энергии 1-10 Дж/см2. Кроме этого, могут предъявляться также определенные требования к однородности плотности потока по сечению пучка, высокой долговечности установки, хорошей повторяемости параметров пучка (обычно нестабильность не должна превышать 10—20%).

С точки зрения вышеперечисленных критериев весьма перспективным представляется использование источников ионов со взрывной ионной эмиссией [4]. В этом случае при энергии ионов Е ~ 100 - 300 кэВ и радиусе пучка г ~ Ю-2 см удается получить плотность выделившейся энергии в интервале Qo ~ 10~2 - 102 Дж/см2, что отвечает плотности выделяемой мощности Ро ~ Ю5 - 108 Вт/см2 при длительности импульса Ai ~ 10~8 -К)"7 с.

Физические модели взаимодействия пучков ионов с веществом

Картина взаимодействия мощных ионных пучков с мишенью весьма сложна. Она содержит большое разнообразие самых различных физических явлений: от тепловых и механических, до газодинамических и плазменных. Хотя многие особенности этих процессов уже ясны, в настоящее время общая физическая модель, описывающая в широком диапазоне параметров такие взаимодействия, пока отсутствует. Характерная черта воздействия мощного пучка ионов на металл состоит в том, что зона структурных изменений в этом случае более чем на два порядка превосходит длину свободного пробега ионов (эффект дальнодействия). Например, дислокационные структуры, образующиеся при ионном облучении, могут располагаться на расстояниях до 100 мкм от поверхности при средней длине свободного пробега ионов менее 1 мкм. Удовлетворительного объяснения этого явления в настоящее время не существует. В частности, гипотеза, согласно которой дефекты образуются вблизи поверхности, а затем мигрируют на большие глубины [5, 6], противоречит общепринятым представлениям о том, что дефекты не могут преодолевать границы зерен вещества. Учет механизмов образования в облучаемых материалах волн сжатия, переходящих в процессе нелинейной эволюции в ударные, представляется весьма привлекательным. При этом дислокации могут образовываться за счет энергии упругих волн, возбуждаемых в образце [7, 8, 9]. Однако и в этом случае имеются экспериментальные данные, которые, по-видимому, не укладываются в рамки такого описания. Мы имеем в виду образование дефектов на расстоянии 100 мкм от источника размером 10 нм, поскольку в этом случае расходящаяся сферическая акустическая волна должна быстро терять плотность энергии, необходимую для образования дефектов [10].

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

При низких интенсивностях пучка, когда средняя тепловая энергия е ~ кТ молекул среды существенно меньше теплоты испарения АН, вылет вещества из мишени может осуществляться только за счет высокоэнергетических хвостов распределений - реализуется тепловой активационный механизм испарения. Скорость распространения фронта испарения V в таком режиме описывается формулой Аррениуса [11]:

V ~ - ехр{-АН/кТ). т

Здесь параметры а - постоянная решетки и г - период собственных колебаний молекул, связанные со скоростью звука с соотношением а/т ~ с. При более высоких мощностях пучка, когда испарение еще не носит характера взрыва, все еще можно говорить о существовании у образца поверхности, разделяющей паровую и конденсированную фазы. Очевидно, что скорость перемещения этой поверхности детерминирована законом сохранения энергии, а учет испарения в уравнении теплопроводности может быть произведен путем наложения соответствующих условий на границе

раздела конденсированной и паровой фаз. В простейшей форме граничное условие на фронте испарения имеет вид [12] дТ

K~ + F~LRVf = О, где К - коэффициент теплопроводности конденсированной среды (теплопроводностью паровой фазы здесь пренебрегаем), F - поглощенная в волне испарения плотность энергии падающего пучка, L - удельная теплота плавления, R - плотность конденсата, Т - его температура, V/ - скорость перемещения фронта испарения. Если Т превышает температуру плавления, то испарение происходит не с твердой, а с жидкой (кипящей) поверхности при давлении насыщенного пара, причем температура пара и вещества с разных сторон поверхности совпадают и равны температуре кипения. Расчет температурных полей в этом случае может быть выполнен с использованием уравнения теплопроводности в рамках задачи Стефана, учитывающей скорость движения границы раздела твердой и жидкой фаз вещества на основе уравнения теплового баланса [13]. При этом температура на границе раздела, как и в случае границы жидкость-пар, меняется непрерывно, а ее производная претерпевает скачок, описываемый соотношением теплового баланса: дт(хь-о,г) дт{хь + о, t)

----~-= QKVb, ох ох где Q - удельная теплота фазового перехода, Kl,Ks ~ теплопроводности в жидкой и твердой фазах соответственно, Vj, - скорость движения границы: Vb = dxb/dt. Динамика фазового перехода типа плавления рассматривается подробно в разделе 1.2.

При дальнейшем повышении плотности энергии падающего пучка механизм равновесного испарения нарушается. В настоящее время описание процессов испарения такого типа производится в рамках модели с кнудсе-новским слоем [И]. Критерием перехода от модели с кипящей поверхностью к модели с кнудсеновским слоем является, очевидно, условие

Р > Р

1 i сг> где Р - давление вблизи поверхности, Pcr ~ критическое давление перехода жидкость-пар. Расчеты, с использованием уравнения Больцмана показывают, что в кнудсеновском слое пара, непосредственно прилегающем к поверхности конденсата и имеющем толщину порядка нескольких длин свободного пробега, устанавливается температура Тк несколько меньшая, чем температура Т жидкой фазы:

Тк = 0,65Т и плотность

РК = Pv, где ру - плотность насыщенного пара. Кроме этого, на границе раздела следует задать условия, отвечающие законам сохранения энергии, импульса и массы.

Даже этот метод описания испарения может стать совсем неадекватным задаче в том случае, когда средняя энергия молекул е в приповерхностной области намного превосходит теплоту испарения. В этих условиях резкой границы между конденсированным и испарившимся веществом не существует, так что уравнения среды в переходной области следует решать совместно с уравнениями состояния. Подобные граничные условия изучались в задачах инерционного ядерного синтеза (см., например, [14]), когда генерируются интенсивные ударные волны, т.е. фактически происходит взрыв поверхностного слоя вещества. Критерии перехода от модели двухфазной границы к этому режиму испарения сформулированы в работе [15].

Фазовые переходы в рамках задачи Стефана

Как известно, моделирование фазовых переходов типа плавления - затвердевания! испарения - конденсации приводит к задаче Стефана. Соответствующие математические модели характеризуются наличием подвижной, заранее неизвестной границы S фазового перехода (см., например, [23]). Формулировка задачи Стефана впервые была дана Г. Ламе и Б.П. Клапейроном в 1831 году, для частного случая сосуществования жидкой и кри-сталической фаз при одинаковой температуре [24]. В 1889 Стефан опубликовал четыре работы, в которых была сформулирована эта проблема в общей форме; и которая, фактически, принята и в настоящее время [25]. В соответствии с этим подходом для межфазовой границы принимается следующее условие

Г{х3 + О, t) дТ{х3 - 0; t) soV-дх--q-дх-= LPsoiVs• (v.3)

Здесь Vs = d^s/dt - скорость движения границы: S, Ksoi и Кцд - коэффициенты теплопроводности материала для твердой и жидкой фаз, L и psoi - удельная теплота плавления (или энтальпия фазового,перехода) и плотность соответственно. Условие (v.3) имеет ясную физическую интерпретацию. Действительно, согласно закону Фурье, поток тепла j пропорционален градиенту температуры, j =- —К grad Т.

Поэтому, левая сторона (v.3), представляющая собой разность потоков тепла поступающего на S и выходящего из нее, описывает тепло, поглощенное единицей площади границы за единицу времени. Выражение в правой, стороне (v.3), очевидно, представляет собой тепло, связанное с плавлением или затвердеванием материала, который за единицу времени пересекла единица площади перемещающейся границы раздела S.

Полная математическая формулировка задачи Стефана включает, кроме (v.3), условие, учитывающее тот факт, что фазовый переход происходит при постоянной температуре,

T|s = TV • (v.4) где Т* - температура плавления, а также - закон сохранения энергии: дТ

PC = -divj + q(x,t).

Здесь g(x, t) описыват мощность внешнего источника тепла, а С - коэффициент теплоемкости. В основополагающей работе Стефана рассматривался частный случай.^(х, t). = 0, так что вся динамика теплопереноса полностью определялась лишь температурным градиентом в среде.

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

Обычно соотношения (v.3) - (v.4) используются в численных алгоритмах в явном виде [26]. В 1953 году А.Н. Тихоновом и А. А. Самарским был предложен другой подход [27], согласно которому условия (v.3) - (v.4) включаются непосредственно в уравнение сохранения энергии. В результате получаем следующую обобщенную формулировку задачи Стефана: рС + L 5{Т - Г*)) ^ + v grad т) = div(K grad Т) + д(х,£), (v.5) где L 5(T — T*) dT/dt дополнительный вклад в теплоемкость тепла, израсходованного на фазовый переход, v grad Т учитывает возможное температурное изменение, обусловленное конвекцией (в задачах, рассматриваемых в данной главе, не возникает необходимости учитывать этот член). Таким образом, главная идея этого подхода сводится к предложению учета удельной теплоты плавления L в качестве дополнительной компоненты теплоемкости рС, которая, однако, дает вклад только в точке фазового перехода, когда Т = Т\

В дальнейшем Самарский и его ученики превратили эту идею в эффективные численные алгоритмы (см., например [28, 29]). Но даже в этом подходе уравнение (v.5) обычно считается только следствием условия (v.3).

В разделе 1.2 будет показано, что на самом деле условие (v.5) дает правильное описание фазовых переходов даже в том случае, когда соотношения (v.3) и (v.4) не применимы.

Процессы трекообразования в раках модели температурного пика (МТП)

Латентные треки ускоренных ионов в высокотемпературных сверхпроводниках способны играть роль центров пиннинга магнитных потоков, существенно увеличивая при этом плотность критического тока [33, 34, 35]. Не смотря на большую практическую важность этих процессов, удовлетворительной теории формирования треков в высокотемпературных сверхпроводниках пока не существует. Теоретически возможны несколько различных механизмов, среди которых модель температурного пика (термопика) [36], модель ионного пика [37], а также другие более изощренные модели (см., например, [38]). Модель ионного пика объясняет формирование трека созданием положительного ионного облака вокруг траектории иона, которое "взрывается" из-за электростатического отталкивания - происходит так называемый кулоновский взрыв. Согласно модели термопика, материал плавится, в пределах цилиндра вокруг траектории ускоренного иона, если температура решетки превышает температуру плавления. Последующий процесс быстрого остывания приводит к формированию аморфной области, которая может наблюдаться методами электронной микроскопии высокого разрешения1.

Первая попытка разработки модели термопика (МТП) для описания процесса формирования треков в высокотемпературных сверхпроводниках была предпринята в [39]. Физическими причинами, которые привели авторов к этому, послужили следующие факты:

• нанодифракция в УВагСизС^-я.показывает, что материал в области треков находится в аморфном состоянии; этом смысле укоренившиеся термины "латентный" или "скрытый" по отношению к рассматриваемым в диссертации трекам не вполне корректны.

• электронная микроскопия показала наличие искажения решетки вокруг треков, соответствующее расширению материала внутри треков [39, 40].

Хотя оба эти факта действительно могут быть интерпретированы естественным образом как последствия процессов плавления материала (сопровождающегося его расширением) и последующего затвердевания, вполне возможны также и другие объяснения. Поэтому исследования с использованием почти той же самой экспериментальной-техники привели авторов работы [34] к другому заключению о том, что механизм формирования трека в УВазСизО^ основан на процессах ионизации [37].

Первое описание формирования треков в высокотемпературных сверхпроводниках с помощью МТП, предпринятое в [39], не учитывало теплоту плавления и поэтому предсказывало существенно большие размеры треков, чем экспериментально измеренные. Чтобы каким-то образом интерпретировать это различие, в [39] была предложена интересная гипотеза "эпитаксиального восстановления", согласно которой не вся расплавленная область переходит в аморфную фазу, а только ее внутренняя часть, в то время как внешняя часть "восстанавливается", т.е. возвращается в прежнее кристаллическое состояние.

В работе [41] для объяснения эволюции размеров треков с изменением энергии налетающего иона для облучаемых сверхпроводников УВагСизОу-^ и Bi2Sr2CaCu208 был предложен феноменологический подход, также основанный на концепции термопика. Хотя эта модель успешно объясняла ряд наблюдаемых особенностей трекообразования, она содержала параметры, не зависящие от физических свойств материалов и, поэтому, ее можно рассматривать только как полезное предварительное исследование проблемы.

Более детальная модель процесса формирования треков в УВагСизОу-^, основанная на системе связанных уравнений для температур электронов и атомов, была предложена в работе [42] по аналогии с моделью термопика [31], развитой в Казне (Франция) для описания треков в аморфных металлах и полупроводниках. В качестве свободного параметра в этой версии МТП принята средняя длина свободного пробега электрона относительно атомов решетки, Л =• \/Der. Здесь De - температуропроводность возбужденных электронов вблизи траектории иона, которая предполагается постоянной (для данного материала) и приблизительно равной 1 - 2 см2/с [43]. Параметр г - время электронно-атомной релаксации - может быть оценен с использованием результатов фемтосекундых лазерных экспериментах [44, 45]. Другие количественные параметры, используемые в модели, это известные теплофизические характеристики облучаемого материала, такие как теплопроводность электронов и атомов Ке и К^ теплоемкость Се и Сг, плотность р твердой и жидкой фаз, температура плавления Тт, удельная теплота плавления Qf.

Значение параметра А ~ 18 нм, найденное в [42] для УВагСизОу-^, оказалось близким к соответствующей величине, полученной для аморфных металлов и полупроводников, а время электронно-атомной релаксации г находилось в хорошем согласии с лазерными фемтосекундными экспериментами. Однако, простые аналитические оценки выполненные в [42], показали, что экспериментально наблюдаемую зависимость радиусов треков от выделившейся при прохождении иона энергии можно объяснить только в том случае, если предположить приблизительно линейную зависимость т от Те (такая зависимость следует, в частности, из теории Аллена [46]). В этом отношении описание формирования трека в УВагСизС^-я оказалось существенно отличным от предложенного в Казне, где величина параметра г считалась не зависящей от температуры. В диссертации т(Те) зависимость учтена в явном виде, т.е. функция т(Те) включена в систему уравнений, описывающих формирование трека.

Охлопывание кавитационного пузырька в D-ацетоне

С тех пор как в 2002 году в работе [78] было, сообщено о регистрации продуктов термоядерных реакций D(d,p)T и D(d,n)3He на,уровне 105 с-1 при акустической кавитации в СзБеО, эта область исследований стала объектом пристального внимания не только акустиков, но и специалистов по ядерной физике, а также физике плазмы. Последние экспериментальные данные по наблюдению продуктов ядерных реакций при схлопывании кавитационного пузырька в D-ацетоне опубликованы в [77]. Они по-прежнему свидетельствуют в пользу этой новой, чрезвычайно простой возможности осуществления реакций термоядерного синтеза2.

Описание процесса схлопывания пузырька, предложенное Р.И. Нигма-тулиным с соавторами [80], основано на численном решении следующей системы уравнений в форме законов сохранения массы, энергии и импульса для газообразной и жидкой фаз в приближении сферической симметрии коллапсирующей кавитационной полости: I + Я (Риг2) = 0, + £| (ри^).+% = о, ч где р,и,р,Т,е и А - плотность массы,-скорость, давление, температура, плотность полной энергии и коэффициент теплопроводности среды.

В решаемую систему входят также уравнения состояния и условие Герца-Кнудсена-Ленгмюра [79] на границе фазового раздела, учитывающее коэффициент конденсации а. Для ацетона а — 1. Это означает, что доля молекул пара, которые при ударе о поверхность раздела фаз отражаются обратно, пренебрежимо мала. По этой причине общая масса пара при сжатии пузырька быстро уменьшается. Кроме того, в случае поглощения молекулы границе раздела передается импульс, который приблизительно в два раза меньше, чем импульс, передаваемый при отражении. Все это сильно снижает сопротивление паров схлопыванию пузырька и, по-видимому,

7 2Результаты, приведенные в [?], имели широкий научный резонанс и стимулировали весьма оживленную дискуссию (см. обзор [84].) является одной из основных причин создания высоких температур и давлений при акустической кавитации в D-ацетоне.

При сверхвысоких давлениях в [80] используются уравнение состояния Ми-Грюнайзена [81], потенциальные функции Борна-Майера [81], которые описывают сверхсжатые жидкости, и современные экспериментальные данные об ударной адиабате в ацетоне [82]. Учитываются диссоциация и ионизация молекул в момент схлопывания пузырька, и соответствующие им потери энергии, а также влияние молекулярной, электронной и ионной теплопроводности.

Численное решение полученной системы уравнений было найдено в [80] с помощью метода Годунова [83]. В результе было установлено, что значения температуры ионов Т и плотности р в пузырьке на стадии схлопывания достигают очень больших значений: Т ~ 107 -г Ю8 К, р ~ 100 г/см3 на расстояниях г ~ 5 -г 25 нм от центра пузырька в течение интервалов времени At ~ Ю-13 -f- Ю-12 с [78].

Выход термоядерных нейтронов в расчете на один пузырек может быть рассчитан по формуле где П£> - объемная концентрация ядер дейтерия; < va(y) > - эффективная скорость реакций D — D-слияния. Он оказывается в хорошем согласии с экспериментом [84].

В данной модели, в соответствии с гипотезой Р.И. Нигматулина, для диссоциированной фазы берутся те же уравнения состояния, что и для недиссоциированного ацетона ("замороженная" ударная адиабата). Хотя подобное предположение и противоречит экспериментальной ударной адиабате [82], оно может быть оправдано быстротечностью завершающей стадии процесса схлопывания. Поскольку характерные времена протекания процесса в данном случае существенно короче, чем в экспериментах, описанных в [82], можно допустить, что существенная перестройка пространственных положений фрагментов, по сравнению с их расположением в молекулах; не успевает произойти. Таким образом,, физический смысл этой гипотезы сводится к приближенному учету термодинамической неравновесности процесса.

Кроме этой гипотезы, в теоретическом описании [80] используются два правдоподобных, хотя и не доказанных строго, предположения о том, что в многопузырьковой системе в акустическом поле во время схлопывания формируется дополнительный положительный импульс давления в несколько десятков атмосфер и о том, что электроны плазмы не успевают нагреться за время схлопывания 10~12 с). Последнее предположение можно назвать гипотезой перегрева ионной компоненты неравновесной плазмы.

Диссертационная работа устроена следующим образом:

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

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

Третья глава посвящена проверке гипотезы перегрева, а также связанному с этим развитию модели в направлении явного учета процессов охлаждения ядер электронами в сверхплотной сильно неравновесной плазме, образующейся при схлопывании кавитационного пузырька в D-ацетоне. Поскольку в данной модели предполагается, что электроны характеризуются определенной температурой Те, расчеты энергии электронов Ее> фактически сводятся к расчетам их теплоемкости Се(Те). Информация о параметре электрон-ионной связи д, как и в процессе трекообразования, эквивалента знанию времени т электрон-ионной релаксации. Основные проблемы, связанные с предлагаемым в диссертации обобщением модели [82] связаны именно с необходимостью расчета параметров g и Се{Те), так как информация об остальных параметрах уже содержится в той или иной форме в работах [82] и [51].

Для исследований процессов термоупругой релаксации описаных в первой главе диссертации и для изучения процессов трекообразования описаных во второй главе, разработан комплекс программ на языке FORTRAN.

Для расчетовов теплоемкости электронов в дейтерированном ацетоне в процессах схлопывания кавитационного пузырька разработан комплекс программ в среде SCILAB 3.0-RC1 (May 17, 2004), которая. является бесплатно распространяем аналогом MATLAB. Тексты программ находятся в приложении.

Заключение диссертация на тему "Математическое моделирование процессов высокотемпературной релаксации в электронно-атомных системах"

3.1.5 Выводы

Оценки с использованием аналитической формулы (3.4), полученной Спитцером [92], показывают, что электронная компонента плазмы действительно может оставаться относительно холодной на протяжении заметной части времени существования сверхплотного состояния вещества в схлопы-вающемся кавитационном пузырьке в CsDgO. Расчеты, выполненные на основе такого модельного допущения, предсказывают в этом случае возможность осуществления термоядерных реакций D-D слияния со скоростью порядка одного события на одно схлопывание пузырька при температуре D-ацетона около 273 К [78, 84]. Оптимизма относительного этого новой возможности реализации термоядерного синтеза прибавляют также два дополнительных обстоятельства. Во-первых, согласно [95], более тщательное численное моделирование показывает, что оценки, выполненные с помощью аналитических формул типа (3.4), использовашихся в данной работе, дают некоторое занижение величины времени релаксации. Во-вторых, расчеты в рамках гидродинамической модели [78, 84], говорят о том, что дальнейшее понижение температуры D-ацетона приводит к уменьшению времени существования термоядерной плазмы почти до значения Ю-13 с, в то время как выход нейтронов при этом увеличивается.

3.2 Теплоемкость электронов в плазме, образующейся при охлопывании кавитационного пузырька в D-ацетоне

3.2.1 Фомулировка модели и результаты расчета

Предполагая локальную компенсацию электрического заряда электронов и ионов, среднее число электронов с энергией Е будем описывать распределением Ферми для невзаимодействующих частиц:

W = е(Е-»(т))/квт + i ' (3-9) где Т - температура электронов, К в - постоянная Больцмана, /л(Т) - химический потенциал электронов, находящихся в электрическом поле ионов. Поскольку для ацетона температура диссоциации Т^ существенно меньше температуры ионизации Ti (Td ~ 3 • 103 К, Т; ~ 1,2 • 105 К), будем считать, что ионизация происходит позже диссоциации и, следовательно, энергии связи электронов совпадают с потенциалами ионизации отдельных атомов [96].

Из того, что при Т = 0 все электроны находятся в связанном состоянии следует, что значение химического потенциала ц при нулевой температуре равно наибольшей энергии связи электрона в атомах С, D и О, т.е. 871,1 эВ. Выделяя в выражении для химического потенциала эту величину отдельным слагаемым,

КТ)=М0) + Дм(Т), выражению (1) нетрудно придать вид е-(а+А»)/квт + i ' е(ек-л цуквТ + х ' соответственно для связанных и свободных состояний. Здесь

Ei = fi(0) - Ei энергия связи, отвечающая г—му электронному уровню, р1

Sk —

2 т энергия свободного электрона. Полученные формулы позволяют трактовать химический потенциал как своего рода потенциальную яму, при выходе из которой электрон становится свободным (затратив на выход энергию, равную энергии связи). Зависящая от температуры поправка к постоянной части химического потенциала А/л(Т) находится из условия сохранения полного числа электронов г к где суммирование выполняется с учетом кратности вхождения атома данного сорта в молекулу СзБбО. Заменяя обычным образом [93] сумму по состояниям свободных электронов интегралом (здесь v - объем, приходящийся на одну молекулу вещества) л/2г>т3/2 v^de- y/2v(mKBTf'2 J.zdz

TT2/i3 J о e(e-An)/KBT+1 JQ ez~Av/KBT + 1 ' получим уравнение для определения величины поправки Afi(T)

Ущ(\- , Л * Т-1 = const Т3/2 Г /ldl-, (3.11)

Z-J г у е-{£г+Ац)/КвТ + IJ JQ ez-An/KBT + I' ^ > где щ - кратность вхождения атомов в формулу вещества.

На рис. 3.2 показана.расчетная зависимость А/л(Т) для СзБбО. Тот факт, что поправка Ац(Т) имеет отрицательное значение объясняется большим статистическим весом свободных состояний, приводящим к понижению температуры начала эффективной ионизации £i/К в, которую можно было бы ожидать без учета этого обстоятельтва. На рис. 3.3 показано среднее число электронов,

NU) = 1 - (Nitb), 82

Т [106 К]

Рис. 3.2. Расчетная зависимость Ац{Т) для СзБ60 при сжатиях р/р0 = 1, 10, 50 и 100 (кривые 1, 2, 3 и 4 соответственно). 1

0.8 ^ 0.6 0.4 0.2

0 0.5 1 1.5 2 2.5

Т [106 К]

Рис. 3.3. Среднее число электронов, выбиваемых с уровней £i — 11.3, 138.1, 490, 871.1 (кривые 1, 2, 3 и 4 соответственно).

Рис. 3.4. Температура 90% ионизации различных электронных уровней для степени сжатия р/ро — 1 и 100.

Рис. 3.5. Температура 90% ионизации в зависимости от сжатия р/ро для самого глубокого электронного уровня £i — 871.1 эВ. освобождающихся с i—го уровня в зависимости от температуры. Видно, что представленная здесь картина последовательного возбуждения все более глубоких уровней вполне аналогична той, которая имеет место при ударной ионизации или в равновесной плазме. Рисунки 3.4 и 3.5 иллюстрируют некоторые важные с физической точки зрения характеристики плазмы, образующейся при схлопывании кавитационного пузырька в D-ацетоне. Так, на рис. 3.4 показана зависимость температуры почти полной ионизации электронных уровней от их энергии связи. Рис.3.5 дает информацию о зависимости температуры почти полной ионизации всех атомов плазмы от степени сжатия.

Теплоемкость электронов, отнесенная к одной молекуле вещества, может быть найдена прямым дифференцированием полной энергии свободных электронов: где суммирование выполняется по всем связанным состояния электронов. На рис. 3.6 представлена рассчитанная таким образом теплоемкость электронов в области сильной ионизации при сжатиях р/ро =1,10 и 100, переведенная в боле'е удобные для практического использования единицы. Видно, что она не описывается линейной зависимостью С = 7Т, справедливой для вырожденного электронного газа в широкой потенциальной яме с плоским дном при температурах меньше или порядка температуры Ферми, когда постоянная (Зоммерфельда) 7 может быть рассчитана по формуле где N - полное число электронов в молекуле [93].

3.2.2 Интерполяционные формулы

Для моделирования процессов схлопывания кавитационого пузырька важно иметь простые, но в тоже время, достаточно точные интерполяци

6 4 2 О О

0.1

0.2

0 0.2 0.4 0.6 0.8

1.2 1.4 1.6 1.8 2

Т [10е К]

Рис. 3.6. Теплоемкость электронов для молекулы D-ацетона при разных сжатиях р/ро, рассчитанная в соответствии с рассматриваемой моделью. онные формулы расчета зависимости Су(Те, р/ро). Их нахождение представляет собой отдельную вычислительную задачу, которой посвящен настоящий раздел. С этой целью в табл. 3.1 приведены соответствующие рис. 3.6 значения теплоемкости Су для некоторых характерных температур электронов. Точки в табл. 3.1 выбирались следующим образом. С помощью программы 1, помещенной в приложении 2, вычислялись значения Су для Те (:Tmin,Tmax) К с точностью Ю-4, где Tmin =9000, 9500 и 10000 К, Ттах = 2 • 106 К для р/ро =1, 10 и 100. Таким образом, для принятого разбиения

Cvi = Cy(Ti), Ti = Tmin т • г, % — 0,., п, т = (Ттах — Tmin)/n, где п - число делений рассматриваемого интервала температур. Потом выбирались предварительные значения элементов таблицы (Т^Сук), к £ (0, п) в точках, для которых первая или вторая производная Су равны нулю, а также некоторые дополнительные промежуточные точки, в которых значения теплоемкости составляли определенный процент относи

1 1 1

1 Ш\АлААа^ Р/Р0=1

0 0.5 1 1.5 2

IML Л/VA/W^ /V\ p/pQ=10 0 0.5 1 1.5 i i i 2 p/pn=100

0.3 о

UJUVVY.А /\

О 0.5 1 1.5 2

T [10е К]

Рис. 3.7. Относительная ошибка интерполяции кубическим сплайном теплоемкости при промежуточных температурах, получающаяся при использовании табличных значений теплоемкости. тельно скачка между экстремальными значениями на соответствующих температурных интервалах. Таким образом, число табулируемых точек для каждого табличного значения степени сжатия p/pQ было фиксированным. На границах интервалов температур в точках Tmin и Ттах вычислялись приближенные значения производной теплоемкости по формулам Cy{Tmin) = (CVi - Cvo)/t и С'у(Ттах) = (CVn - Cvn-i)/T. Эти значения приведены в табл. 3.2.

Далее по выбранным точкам строился кубический сплайн для теплоемкости с вышеуказанными значениями производной на границах, с целью сравнения полученных таким образом величин Су к с точными. Это позволило определить характерную для табл. 3.1 относительную ошибку интерполяции теплоемкости при промежуточных, не охватываемых таблицей, значениях температур электронов. Далее предварительные элементы

Заключение

В дисертационной работе проведены исследования процессов термоупругой релаксации и процессов трекообразования в материаллах облучаемых ускоренными ионными пучками. Разработаны модели, алгоритмы и программы для математического моделирования этих процессов и для вычисления неизвестных параметров моделей. Основные результаты работы можно сформулировать в следующих пунктах:

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

2. Разработан комплекс математических моделей и программ для расчета процессов тепловой релаксации в образцах, подвергаемых облучению ускоренными ионами. С его помощью были выполнены исследования взаимодействий ионов с веществом и получены следующие результаты: а) Изучен процесс формирования и релаксации термоупругих напряжений в материалах, испытывающих интенсивную ионную бомбардировку. Сформулирована модель пространственно-временной динамики энерговыделения в тонких пленках, отвечающая имеющимся источникам со взрывной ионной эмиссией. Установлено, что максимальные механические напряжения, возникающие вследствие теплового расширения вещества, не превышают ни предел текучести железа, ни его временное сопротивление разрыву и таким образом не являются ответственными за модификацию образца. б) В результате численных экспериментов и аналитических оценок получен вывод о невозможности постановки традиционной задачи Стефана для описания движения границы раздела твердой и жидкой фаз в том случае, когда в уравнении модели присутствует источник энерговыделения. в) В рамках концепции температурного пика рассмотрены процессы трекообразования в высокотемпературных сверхпроводниках. С учетом имеющихся теоретических и экспериментальных данных разработана модель процесса энеговыделения, описывающая его пространственно-временную динамику. Уточнена модель электронной подсистемы с использованием теории Аллена и результатов фемтосекундных лазерных экспериментов. г) Продемонстрирована применимость концепции температурного пика для описания процессов трекообразования в иттриевых сверхпроводниках. Установлена бифуркационная зависимость решений от параметра температуропроводности электронов, связанная с процессом электронной закалки.

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

Благодарности

Выражаю огромную благодарность своим научным руководителям Б.Ф. Костенко и Э.А. Айряну за постановку задач, поддержку и постоянное внимание. Выражаю глубокую признательность профессорам Е.П. Жидкову и И.В. Пузынину за проявленный интерес к работе и полезные дискуссии. Особо хочу поблагодарить О.И. Стрельцову, Я. Бушу, И. Покорны за ценные консультации. Так же хочу поблагодарить соавторов и коллег по работе: И.В. Амирханова, И.Н. Гончарова, Т.А. Стриж, Д.В. Подгайного и Т. JI. Бояджиева. Искренне благодарен своей супруге Юдите за неоценимую помощь и поддержку.

Отдельно хочу поблагодарить дирекцию Лаборатории информационных технологий Объединенного института ядерных исследований за предоставленные хорошие условия для работы.

Библиография Прибиш, Ян, диссертация по теме Математическое моделирование, численные методы и комплексы программ

1. В.И. Бойко, А.Н. Валиев, А.Д. Погребняк, Модификация металлических материалов импульсными мощными пучками частиц, УФН, 1999, т.169(11), с.1243.

2. А.Н. Диденко и др., Поверхность, 1985, вып. 1, с.150.

3. А.N. Didenko et а/., VI Int. Conf. on Effect, and Ion Beam, Kobe, Japan, 1986.

4. C.A. Коренев и др. Труды I Межд. симпозиума по пучковым технологиям (ВТ'95), Дубна, 1995.

5. Н.П. Морозов, Д.И. Тетельбаум, Связь аморфизуемости алмазо-подобных полупроводников с их механическими свойствами, ФТП, 1983, т. 17, с. 838.

6. Yu. Martynenko, P. Moscovkin, Dislocation structure formation due to ion bombardment, Rad. Eff. and Defects in Solids, 1994, v. 129, p. 193.

7. V.P. Zhukov, A: V. Ryalenko, Rad. Eff., 1984, No. 1-2, p. 85.

8. Семин Ю.А., Скупое В.Д., Тетельбаум Д.И., Письма в ЖТФ, 1988, т. 14, с. 273.

9. Павлов П. В., Семин Ю. А., Скупое В Д., Тетельбаум Д. Я., Влияние упругих волн, возникающих при ионной бомбардировке, на структурное совершенство полупроводниковых кристаллов, ФТП, 1986, т. 20, с. 503.

10. Ю.В. Мартыненко, Эффекты дальнодействия при ионной имплантации, в сб.: Итоги науки и техники, сер. Пучки заряженных частиц и твердое тело, ВИНИТИ, Москва, 1993, т. 7, с. 82.

11. С.И.Анисимов и др., Действие излучения большой мощности на металлы, М.: Наука, 1970.

12. В.И. Зубов и др., Журн. вычислит, мат. и мат. физ., 1986, т. 26, с. 1740.

13. Айрян Э.А. и др.у Численное моделирование термальной обработки поверхности металла сильноточным ионным пучком, JINR Rapid Comm. 1997, No.686.-97, с. 103

14. Аккерман А.Ф., Бушман А.В., Демидов Б.А. и др., ЖЭТФ, 1985, т.89, с. 852.

15. В.П. Кривобокое, О.В. Пащенко, Г.А. Сапульская, Компьютерное моделирование эрозии и термомеханических процессов в твердом теле, облучаемом мощными наносекундными пучками заряженных частиц, Изв. вузов. Физика, 1993, 12, с.37.

16. И.В. Амирханов и др., Препринт ОИЯИ Р2-98-63; ОИЯИ Р2-98-201, Дубна, 1998.

17. В. Новацкий, Вопросы термоупругости. Москва, АН СССР, 1962.

18. С. А. Коренев, Препринт ОИЯИ Р13-89-615, Дубна, 1989.

19. L.C. Northcliffe, R.F. Schilling, Nuclear Data Tables, 1970, A7, p. 233.

20. J.P. Biersack, L.G. Haggmark, A Monte Carlo computer program for the transport of energetic ions in amorphous targets, Nucl.' Instr. and Meth., 1980, v. 174, p. 257.

21. А.А. Самарский, Теория разностных схем. Москва, Наука, 1977.

22. Н.Н. Калиткин, Численные методы. Москва, Наука, 1978.

23. A.A. Samarskii, P.N. Vabishchevich, Computational Heat Transfer, Mathematical Modelling, V.l, p. 30-33 (John Wiley k Sons, Chichester New-York), 1995.

24. G. Lame, B.P. Clapeiron, Ann. de Chem. et de Phys. XLVII, 1831, pp. 250-256.

25. J. Stefan, Sitzber. Wien. Akad. Mat. naturw. 98, 1889, pp. 473-484, 616634, 965-983, 1418-1442.

26. И.В. Амирханов, E.B. Земляная, И.В. Пузынин, Т.П. Пузынина, И. Сархадов, Численное моделирование фазовых переходов в металах, облучаемых импульсными пучками ионов, Сообщение ОИЯИ, Р11-2001-164.

27. А.Н. Тихонов, А.А. Самарский, Уравнения математической физики. Москва, Наука, 1953.

28. А.А. Самарский, B.D. Moiseenko, Журн. вычислит, мат. и мат. физ., 5, 1965, с. 816.

29. А.А. Самарский, П.Н. Вабищевич, Вычислительная теплопередача. УРСС, Москва, 2003.

30. А. А. Самарский, А.В. Гулин, Численные методы математической физики. Научный мир, Москва, 2000.

31. Г. И. Марчук, Методы вычислительной математики. Москва, Наука, 1977

32. Р.А. Влейхер, В.П. Кривобокое, С.В. Пащенко, Тепломассоперенос в твердом теле под действием мощных пучков заряженных частиц. Новосибирск, Наука, 1999.

33. M. Toulemonde, С. Dujour, E. Paumier, Phys. Rev. B46, pp. 14362 -14369 (1992-II).

34. A.M. Мейрманов, Задача Стефана. Новосибирск, Наука, 1986.

35. М. Konczykowski, F. Rullier-Albenque, E.R. Yacoby, A. Shaulov, Y. Yeshurun, P. Lejay, Phys. Rev. B44, 1991, p. 7167.

36. V. Hardy, D. Groult, M. Hervieu, J. Provost, B. Raveau, S. Bouffard, Nucl. Instr. and Meth. B54, 1991, p. 472.

37. L. Civale, A.D. Marwick, T.K.Worthington, M.A.Kirk, J.R. Thompson, L. Krusin-Elbaum, Y. Sun, J.R. Clem, F. Holtzberg, Phys. Rev. Lett. 67, 1991, p. 648.

38. I.M. Lifshits, M.I. Kaganov, L. V. Tanatarov, J. Nucl. Energy A 12,1960, p. 69.

39. R.L. Fleischer, P.B. Price, R.M. Walker, Nuclear Tracks in Solids, University of California Press, Berkeley, 1979.

40. E. Dartyge, P. Sigmund, Phys. Rev. В 32, 1985, p. 5429.

41. Yimei Zhu, Z.X. Cai, R.C. Budhani, M. Suenaga, D.O. Welch, Phys. Rev. В 48, 1993, p. 6436.

42. Yimei Zhu, Z.X. Cai, D.O. Welch, Phil. Mag. A 73, 1996, p. 1.

43. G. Szenes, Phys. Rev., B54, 1996, p. 12458.

44. I.N. Goncharov, B.F. Kostenko, V.P.Philinova, Phys. Lett., A 288, 2001, p. 111.

45. A. Meftah, F. Brisard, J.M. Constantini, E. Dooryhee, M. Hage-Ali, M. Hervieu, J.P. Stoquert, F. Studer, M. Toulemonde, Phys. Rev., B49, 1994, p. 12457.

46. S.D. Brorson, A. Kazeroonian, D.W. Face, Т.К. Cheng, G.L. Doll, M.S. Dresselhaus, G. Dresselhaus, E.P. Ippen, T. Venkatesan, X.D. Wu, A. Inam, Solid State Comm., 74, 1990, p. 1305.

47. I.I. Vengrus, A.L. Dobryakov, S.A. Kovalenko, V.S. Letokhov, U.E. Lozovik, G. Marovski, U.A. Matveets, V.M. Farztdinov, N.R. Ernsting, Pis'ma v JETF, 62, 1995, p. 739.

48. P.B. Allen, X. Du, L. Mihaly, L. Forro, Phys. Rev., B49, 1994, p. 9073.

49. A.E. Volkov, V.A. Borodin, Nucl. Instr. and Meth. in Phys. Research B193, 2002, pp. 381-390.

50. L.T. Chadderton, I.M. Torrens, Fission Damage in Crystals. Methuen, London, 1969.

51. W.C. Moss, D.A. Young, J.A. Harte, J.L. Levatin, B.F. Rozsnyai, G.B.Zimmerman, I.H. Zimmerman, Phys. Rev., E49, 1999, p. 2986.

52. M.P.R. Waligorski, R.N. Hamm, R. Katz, Nucl. Tracks Radiat. Meas., 1, 1986. p. 309.

53. W.H. Barkas Nuclear Research Emulsions I. Techniques and Theory. Academic Press, New York, 1963.

54. Л.Д. Ландау, И.М. Лифшиц. Теория поля. Москва, Наука, 1988.

55. J.F. Ziegler, SRIM 2003, version 2003.26, www.srim org.

56. V.S. Barashenkov, Russ. Chem. of High Energies, 28, 1994, p. 229.

57. A. Meftah, F. Brisard, J.M. Gonstantini, M.Hage-Ali, J.P. Stoquert, F. Studer, M. Toulemonde, Phys. Rev., B48, 1993, p. 920.

58. G. Sciwietz, G. Xiao, E. Luderer, P.L. Grande, Nucl. Instr. and Meth. B164-165, 2000, p. 354.

59. E.A. Ayrjan, A.V. Fedorov, B.F. Kostenko, Part, and Nucl. Lett., 99, 2000, p. 42.

60. H. Krakauer, W.E. Pickett, R.E. Cohen, J.Supercond. 1, 1998, p. 111.

61. M.F. Crommie and A. Zettl, Phys. Rev., B41, 1990, p. 10978.

62. J.L. Cohn et al, Physical Properties of High Temperature Superconductors, III, ed. by D.M. Ginsberg, World Scientific, Singapore, 1992.

63. Yu.V. Martynenko, Yu.N. Yavlinski, Sov. Phys. Dokl., 28, 1983, p. 391; Preprint IAE-4084/11, Moscow, 1985.

64. V.E. Lusternik, V.E. Peletsky, V.S. Bakunov, A.B. Bolotnikov, R\issmn Superconductivity: Physics, Chemistry, Engineering, 3, 1990, p. 2037.

65. H. Teichler, Phys. Rev. B59, 1999, p. 8473.

66. G.M.Mironova, Material Science Forum, 133-136, 1993, p. 847.

67. K. Salama, D.F. Lee, Supercond. Sci. Technol, 7 1994, p. 177;

68. S. Williamson, G. Mourou, J. CM. Li, Phys. Rev. Lett. 52 1984, p. 2364.

69. Y. Idemoto, K. Fueki, Jpn. J. Appl.Phys., 29, 1990, p. 2729.

70. S. Wermbter, L. Tewordt, Physica, C183, 1991, p. 365.

71. S.D. Peacor, J.L. Cohn, C. Uher, Phys. Rev., B43, 1991, 8721.

72. D.R. Atthey, J. Inst. Maths Applies, 13, 1974, p. 353.

73. B.F. Kostenko, J. Pribis, I.V. Puzynin, e-print: math-ph/0302044, 2003, to be published in J. of Comput. Meth. in Sciences and Engineering.

74. F. Faupel, W. Frank, M.-P. Macht, H. Mehrer, V. Naundorf, K. Ratzke, H. R. Schober, S. K. Sharma, H.Teichler, Rev. Mod. Phys., 75, 2003, p. 237.

75. E.P. Donovan, F. Spaepen, D.Turnbull, J.M. Poate, D.C. Jacobson, J. Appl. Phys., 57, 1985, p. 1795.

76. C. Dufour, A. Audouard, F. Beuneu, J. Dural, J.P. Girard, A. Hairie, M. Levalois, E. Paumier, M. Toulemonde, J. Phys.: Condens. Matter, 5, 1993, p. 4573, 1993.

77. Kostenko B.F., Pribis J., Математическое моделирование трекообразования в высокотемпературных сверхпроводниках, Вестник РУДН, сер. Прикладная математика, 2005, т. 4, No 1, с.75-87.

78. Kostenko B.F., Pribis J., Goncharov I.N., Thermal spike model of track formation in УВагСизОг^, Preprint JINR E17-2005-61, Particle and Nuclei Letters, 1(130), Vol. 3, 2006, pp. 31-44.

79. Taleyarkhan R.P., Cho J.S., West C.D., Lahey Jr. R. Т., Nigmatulin R.I., Block R.C., Phys. Rev. E69, 2004, p. 036109.

80. Taleyarkhan R.P., West C.D., Cho J.S., Lahey Jr. R. Т., Nigmatulin R.I., Block R.C., Science 295, 2002, pp. 1868-1873.

81. Schrage R. W., A Theoretical Study of Interphase Mass Transfer. N.Y.: Columbia U.P., 1953.

82. Nigmatulin R.I., Lahey Jr. R.T., Taleyarkhan R.P., West C.D., J. Acoust. Soc. Am 113, 2003, pp. 2205-2206.

83. Зельдович Я.В., Райзер Ю.П., Физика ударных волн и высокотемпературных гидродинамических явлений. Москва, Наука, 1966.

84. Р.Ф Трунин, Л.Ф. Гударенко, М.В. Жерноклетов, Г.В. Симаков., Экспериментальные данные по ударно-волновому сжатию и адиабатическому расширению конденсированных веществ, Саров, РФЯЦ-ВНИИЭФ, 2001.

85. С. К. Годунов, А. В. Забродин, М.Я. Иванов, А.Н. Крайко, Г.П. Прокопов, Численное решение многомерных задач газовой динамики. Москва, Наука, 1976.

86. Беляев В.В., Костенко Б.Ф., Миллер М.Б., Сермягин А.В., Тополъ-ников А. С., Свервысокие температуры и акустическая кавитация, Сообщение ОИЯИ P3-2003-214, Дубна, 2003.

87. Б.Ф. Костенко, Я. Прибиш., Оценка времени охлаждения ядер электронами в сверхплотной сильно неравновесной плазме, Сообщение ОИЯИ Р4-2004-42, Дубна, 2004.

88. Б.Ф. Костенко, Я. Прибиш., Теплоемкость электронов в плазме, образующейся присхлопывании кавитационного пузырька в D-ацетоне, Сообщение ОИЯИ Р11-2004-139, Дубна, 2004.

89. D.J. Flannigan, KS. Suslick, Nature 434, 2005, p. 52.

90. Ландау JI.Д., Лифшиц Е.М., Статистическая механика. Москва, Наука, 1976, с. 186.

91. Robinson В.В, Bernstein L.B., Ann. Phys. (NY), 1962, V. 18, pp. 110169.

92. Гардинер К.В., Стохастические методы в естественных науках. Москва, Мир, 1986.

93. Синельников К.Д., Руткевич Б.Н., Лекции по физике плазмы. Харьков, Изд-во ХГУ, 1964.

94. Спитцер Л., Физика полностью ионизированного газа, Москва, ИЛ, 1957.

95. Л.Д. Ландау, Е.М. Лифшиц, Статистическая физика, 4.1, Москва, Наука, 1976.

96. Н.С. Бахвалов, Н.П. Жидков, Г.М. Кобельков, Численные методы, Москва, Наука, 1987.

97. Клеммоу Ф., Доуэрти Дж., Электродинамика частиц и плазмы. Москва, Мир, 1996, с.434.

98. Справочник физических величин. Под ред. И.К. Кикоина, Москва, Атомиздат, 1976.