УДК 548.4:536.48
АЛГОРИТМИЗАЦИЯ ИМИТАЦИИ ОБРАЗОВАНИЯ ДИСЛОКАЦИОННОЙ ПЕТЛИ ИСТОЧНИКОМ И ПРОЦЕССА ЕЕ ЭВОЛЮЦИИ В ПЛОСКОСТИ КРИСТАЛЛОГРАФИЧЕСКОГО СКОЛЬЖЕНИЯ СО СЛУЧАЙНО РАСПРЕДЕЛЕННЫМИ ДИСКРЕТНЫМИ ПРЕПЯТСТВИЯМИ
Слободской М.И., Матющенко А.В., Голосова Т.Н. (Томск)
Abstract
It is shown that when describing an elementary crystallographic slip a detailed considering of dislocation source action and planar dislocation loop development is necessary instead of using of traditional model of infinitely long straight (or quasi-straight) dislocations. The algorithms for computer simulation of planar dislocation loop emission by Frank-Read source and offurther evolution of the dislocation loop in a random field ofpoint-like obstacles is suggested. The obstacles by spectrum of their “strength ” correspond to reactive and nonreactive forest dislocations.
Широкое применение представлений о дислокациях, сформулированных
Вольтерой в 1907 г., началось после того, как ряд исследователей ( Френкель, Орован,
Поляни, Тейлор и др.) установили связь между дислокациями и пластичностью. После возникновения идеи винтовой дислокации (Бюргерс, 1939), источника дислокаций, способного воспроизводить серию дислокационных петель (Франк и Рид, 1950) и экспериментального наблюдения дислокаций (Хирш, 1956 г.) теория дислокаций и теория пластичности на несколько десятилетий становятся синонимами.
Несмотря на продолжительный бум теоретических и экспериментальных
исследований дислокаций, несомненные успехи в исследовании структуры деформируемых кристаллических тел, значительного продвижения в понимании пластической деформации как целостного явления не произошло.
Одна из причин сложившейся ситуации состоит в том, что многие задачи физики дислокаций решены каждая в своих упрощающих приближениях и с учетом лишь черт, существенных для конкретной предложенной задачи. Существующие знания о дислокационной подсистеме, её эволюции в процессе деформации содержат
многочисленные упущения, неточности и даже некорректности, накапливающиеся в процессе синтеза таких знаний для описания макроскопического поведения кристаллов. Поэтому прямой синтез результатов различных теоретических исследований в единую теорию вряд ли возможен, для построения действительно последовательной и работоспособной теории дислокаций необходима переформулировка большинства уже решенных частных задач.
Так, например, в серии попыток перейти от дислокационных представлений и механизмов непосредственно к макроскопическому поведению кристалла (Тейлор,1934; Мотт,1952; Фридель,1954; Зегер с соавт.,1954-1960; и др.) имеется общее упрощающее приближение, приемлемость которого весьма сомнительна, - приближение
прямолинейных бесконечных дислокаций. В тех случаях, когда речь идет о взаимодействии скользящей дислокации с дислокациями леса и о скользящих
дислокациях, содержащих пороги, допускаются локальные прогибы дислокаций и говорят о квазипрямолинейных дислокациях.
Поскольку в действительности дислокация, не выходящая на поверхность, замкнута, очевидно, что приближение бесконечных прямолинейных дислокаций должно было породить и породило чисто модельные эффекты и проблемы и вместе с тем привело к утрате части информации о дислокационной подсистеме. Такова, например, проблема генераций дислокаций в процессе деформации. В случае представления дислокационной подсистемы бесконечными прямолинейными дислокациями генерацию приходится вводить в теорию как параметр, в то время как при конкретном рассмотрении
распространения кристаллографического сдвига как расширении дислокационной петли деформация и накопление дислокаций оказываются связанными необходимым кристаллогеометрическим соотношением (Инденбом В.Л., Орлов А.Н., 1962).
Приближением бесконечных прямолинейных дислокаций порождено
представление о неразрешимости проблемы unzipping'а - нахождения числа стопоров, пройденных дислокацией атермически, после одной термической активации. Нет даже каких-либо подходов к ее теоретическому решению.
Мы ограничим этот список еще одним следствием приближения дислокационных подсистем прямолинейными бесконечными дислокациями, которое принципиально
меняет всю картину дислокационной пластичности кристалла. Бесконечная
прямолинейная дислокация находится в безразличном равновесии. Ничего не изменится, если поступательно переместить ее на некоторое расстояние Ax. Такие представления в какой-то мере приемлемы для рассмотрения дислокационных подсистем в статическом приближении. Если же учесть, что в действительности дислокации образуют замкнутые дислокационные петли, картина кардинально меняется. Дислокационная петля, если даже предположить, что ее линейное натяжение и радиус кривизны постоянны и она имеет форму, близкую к окружности, может находиться в неустойчивом равновесии лишь при строго определенном напряжении. При малом увеличении радиуса AR петля будет ускоренно расширяться; при малом уменьшении AR - ускоренно сжиматься до аннигиляции.
Таким образом, применительно к дислокационной подсистеме, состоящей из конечных замкнутых петель, статический подход не приемлем. Эволюция дислокационной подсистемы происходит как совокупный результат многочисленных локальных последовательных неустойчивостей, которые сопровождаются движением, размножением и взаимодействием дислокаций в динамическом режиме. В областях локальных потерь устойчивости эволюция дислокации происходит по законам дислокационной микромеханики, адекватно учитывающей кристаллографию, топологию и инерционные свойства дислокаций, сосредоточенные и распределенные силы, действующие на дислокацию, включая силы самодействия дислокаций. Развитие событий в этих областях имеет определенные характерные масштабы и времена и не зависит от длительности деформирующего воздействия.
Макроскопическая деформация складывается из событий, происходящих в областях потери устойчивости. Ее минимальным структурным элементом являются именно такие области, а не дислокации, как в приближении прямолинейных бесконечных дислокаций, характерном для традиционной дислокационной теории пластичности. Но минимальным объектом, который отражает свойство дислокационной подсистемы, является дислокационная петля; именно петля, а не прямолинейная дислокация, которая,
как мы видели, не отражает фундаментальных свойств дислокационной подсистемы - ее способности размножаться в результате многочисленных неустойчивостей. Поэтому первым шагом последовательной теории пластичности к описанию макроскопического пластического поведения кристалла является изучение методами дислокационной микромеханики, включая её статические, динамические и кинематические аспекты, работы дислокационного источника, формирования первой дислокационной петли, испускания серии дислокационных петель.
Однако здесь мы сталкиваемся с принципиальными трудностями. Прежде всего, скользящая дислокация взаимодействует с огромным числом препятствий. Уже число пересечений с дислокациями некомпланарных систем скольжения составляет десятки тысяч, а, кроме того, дислокации взаимодействуют с атомами растворенных элементов, частицами второй фазы и т. д. Разумеется, взаимодействие дислокации со стопором оказывает влияние на ее взаимодействие с другими стопорами; то есть задача описания движения дислокации, осуществляющей кристаллографический сдвиг, является проблемой многих тел, которая относится к числу неинтегрируемых [14]. Кроме того, учёт взаимодействия дислокации с ансамблем случайно расположенных препятствий ставит серьезные теоретико-вероятностные вопросы. Между тем, без описания эволюции дислокационной петли невозможно двигаться дальше к описанию макроскопической пластической деформации кристалла под деформирующими воздействиями.
Приемлемым выходом из этой ситуации представляется описание распространения кристаллографического скольжения методами имитационного моделирования. Существует ещё одна причина принципиальной необходимости использования вычислительного эксперимента при описании кинетики деформационных дефектов : действительные деформационные дефекты структуры, формирующиеся в процессе деформации и осуществляющие ее, существуют в динамике и в сложных взаимодействиях, которые могут быть воспроизведены лишь посредством вычислительного эксперимента на основе математических моделей. Экспериментатор же имеет дело с реликтовыми структурами, сохранившимися в материале после релаксационных процессов динамического возврата.
Методы имитационного моделирования хорошо разработаны для описания движения квазипрямолинейных участков дислокаций по относительно небольшому полю случайно расположенных стопоров (см. обзоры [1,2]). Применение их для описания работы дислокационного источника и эволюции замкнутой планарной дислокационной петли в плоскости кристаллографического скольжения, содержащей препятствия для движущейся дислокации, и для описания макроскопической пластичности предъявило ряд требований к этоим методам. Оказалось [3], что удовлетворительное описание макроскопического пластического поведения кристалла может быть достигнуто уже в рамках трехуровневой модели с использованием единственного мезоскопического структурного уровня - зоны кристаллографического сдвига. Даже в этом случае необходим учет сотен тысяч препятствий [4]. Кроме того, препятствия могут существенно отличаться уровнем напряжений, необходимым для их преодоления скользящей дислокацией, например, реагирующих и нереагирующих дислокаций некомпланарных систем скольжения [13]. Всё выше сказанное ставит в разряд весьма актуальных задачу моделирования процесса образования дислокационной петли источником, ее эволюции в достаточно представительном случайно распределенном поле точечных препятствий с существенно различными прочностями. В процессе алгоритмизации модели появляется круг задач, которые условно можно сгруппировать как вопросы, связанные: а) с эффектом
перехода от представления дислокационной подсистемы бесконечными прямолинейными дислокациями к рассмотрению распространения кристаллографического сдвига как расширения дислокационных петель, сгенерированных источником дислокаций; б) с эффектом “большого” массива препятствий; в) с наличием в плоскости
кристаллографического скольжения препятствий со значительно отличающимися прочностями. Основные трудности в вопросах "а" связаны с неоднозначностью обратных тригонометрических функций, которая приводит к сложностям при идентификации области поиска препятствий в процессе продвижения элемента дислокации и вычислении углов огибания на препятствиях. Требуются критерии, позволяющие отличать верхнюю и нижнюю области на рис.1а, конфигурации "б" и "в" на этом же рисунке.
Рис. 1. Схема возможных конфигураций элемента дислокации
около препятствий
Круг подобных вопросов разрешается с помощью введения определений положительной и отрицательной полуплоскостей относительно вектора СА (рис. 1г) [9] и критерия принадлежности любой точки плоскости, скажем В, указанным полуплоскостям. Нетрудно проверить, что все определяется знаком векторного произведения векторов СА и СВ, а именно:
-1 &В еР-А
sign[CA хСВ ] =
+1 &В ЄР+А 0 <^В єпрямой СА
Присутствие в плоскости кристаллографического скольжения препятствий различной мощности приводит к возможности самопересечений дислокации. Кроме этого, требуется локализовать место замыкания дислокационного сегмента в дислокационную петлю (рис. 1д). Алгоритм замыкания петли в случае препятствий слабой прочности довольно прост и основан на анизотропно-равномерном расширении дислокационного сегмента. Очевидно, что в этом случае замыкание происходит в полосе, расположенной
ниже отрезка МК Поэтому в случае непопадания элемента дислокации в эту область поиск точки замыкания не производится. В противном случае участки МК и КР запоминаются, определяются абсциссы крайнего правого на участке МК и крайнего левого на участке КР стопоров. Пусть это х1 и х2 соответственно. Если | хгх2 I >5, где 5 -константа, определяемая приложенным напряжением, то пересечение участков МК и КР невозможно. В противном случае производится углубленная проверка смежных дуг на пересечение и локализуется место возможного пересечения. В случае пересечения определяются координаты точки пересечения, образовавшаяся петля отшнуровывается. При наличии в плоскости кристаллографического скольжения препятствий с существенно различными прочностями проверку на замыкание петли приходится проводить при попадании дислокационных сегментов в полуплоскость, расположенную ниже прямой, проходящей через точки М и N. Для поиска мест самопересечений дислокации каждый дислокационный сегмент проверяется на возможное пересечение со всеми оставшимися сегментами. Разумеется, соответствующая процедура организована таким образом, что число подобных проверок минимально. На первом этапе также выполняется упрощенная проверка, суть которой состоит в том, что на пересечение проверяются не дислокационные дуги, а соответствующие дислокационные сегменты. При этом активно используется простой критерий пересечения двух отрезков: отрезки ВС и АР
пересекаются тогда и только тогда, когда у векторных произведений [ВСхВА] и [ВСхВР]; [РАхРВ] и [РАхРС] разные знаки [11].
На вопросах, появляющихся вследствие большого массива препятствий, остановимся более подробно.
О р г а н и з а ц и я м а с с и в а п р е п я т с т в и й. Оригинальное решение проблемы “большого массива препятствий” предложено в работе [5] - метод “посыпания” стопоров, сущность которого состоит в том, что препятствия “подбрасываются” по фронту движущейся дислокации. Однако этот метод при всех его несомненных достоинствах не “помнит” массива пройденных стопоров, тогда как в ряде вопросов это необходимо [4,9,11,12]. Кроме того, таким образом сформированное поле препятствий трудно поддается статистическому тестированию на случайность. В настоящей работе число стопоров в плоскости скольжения является параметром модели, ограниченным только временем работы ЭВМ, а сама плоскость кристаллографического скольжения условно разбита на т равных по площади и определенным образом занумерованных блоков. Размер одного блока выбирается исходя из требования минимизации времени работы ЭВМ по обработке дислокационного сегмента и условия, чтобы любой дислокационный сегмент целиком помещался на участке плоскости кристаллографического скольжения не более чем из четырех смежных блоков. Тогда при обработке дислокационного сегмента возможны следующие ситуации (события):
А21-1={для обработки дислокационного сегмента требуется 1
блоков и один из них находится в оперативной памяти
ЭВМ };
А21 = {требуется 1 блоков и ни одного из них нет в памяти},
1=1,2,3,4. Если через Т; обозначить время обработки 1-й ситуации, то
Т21-1=Ь2 [1-Уп+(1-1)Уе], Т21 = 1-Ь2(Уп+Уе), (1)
где Ус - скорость формирования препятствий в блоке, Уп - скорость их просмотра в этом блоке, И - линейный размер блока. Среднее время <Т> обработки одного дислокационного сегмента
Т =1рТ . (2)
¿=1
В (2) р;=Р(А;) - вероятность события А;. Вычисляя р; и минимизируя выражение (2) по И с учетом (1), находим оптимальный линейный размер блока
Иопт = ¿тах (1-Ус/2Уп)/2, (3)
где ётах - максимальная длина дислокационного сегмента, Ус и Уп - упомянутые выше константы конкретной ЭВМ, причем, Ус << Уп .
Заполнение отдельного блока осуществляется отдельной подпрограммой, которая по условно присвоенному этому блоку номеру генерирует число препятствий п в нём в соответствии с распределением Пуассона
Р(п)=(Ь27п!)- ехр(-И2),
формирует координаты препятствий согласно расположению блока в глобальной системе координат по формулам :
Х=АХ|+ ; АХ|=Ь{]-1-1КТ[(]-1)/к]к};
У=АУ|+ ; АУ|=ЬШТ[(]-1)/к]. (4)
Здесь (Х,У), (АXj АУ|) - соответственно глобальные и локальные координаты препятствий ]-го блока, ЮТ(7) - оператор выделения целой части числа 7, £ - равномерно распределенное на интервале [0,1) случайное число. В этой же подпрограмме методом Монте-Карло разыгрывается тип препятствий в соответствии с их прочностями и относительными концентрациями.
Такое поблочное формирование случайного поля стопоров предъявляет особые требования к эффективности датчиков псевдослучайных чисел. В литературе описано несколько “плохих” генераторов случайных чисел [6]. Поэтому в работе, согласно рекомендациям, вытекающим из теории Д.Кнута [7], разработан датчик псевдослучайных равномерно распределенных чисел на интервале [0,1), который, как показали специальные исследования [9], “проходит” все существенные для нашей работы статистические тесты : на равномерность (по критериям х2 , ю2 , Смирнова-Колмогорова), на последовательную корреляцию, проверку подпоследовательностей, проверку серий, проверку на
монотонность, оценку методом Монте-Карло некоторых математических констант [8].
П р е д с т а в л е н и е д и с л о к а ц и и в ЭВМ. Для получения необходимой информации посредством имитации эволюции дислокационной петли в плоскости кристаллографического скольжения, содержащей препятствия, важно, чтобы данные о дислокации хранились в ЭВМ в легкодоступной форме, были в определенном смысле полными, могли бы быть легко модифицированы [10]. Зная координаты точек закрепления дислокации (стопоров, с которыми контактирует дислокация), можно получить всю информацию о дислокации [9]. Поэтому дислокация может быть представлена информацией о точках закрепления, хранящейся в некотором
упорядоченном массиве. Но при таком простейшем представлении дислокации много времени занимает процедура адаптации списка после включения или исключения из него стопора. Поскольку подобная операция в процессе моделирования эволюции одной дислокационной петли встречается сотни тысяч раз, вопрос о представлении дислокации становится особенно актуальным. В предлагаемой методике моделирования эволюции дислокационной петли в плоскости кристаллографического скольжения данные о дислокации хранятся в двух массивах : в одном из них представлены номера и координаты точек закрепления, их прочности и некоторая вспомогательная информация, в другом - номер (присвоенный ему в первом массиве) ближайшего “соседа” на дислокации против и по часовой стрелки. Можно показать, что адаптация такого списка производится с минимально возможным числом перевычислений.
Организация активационного процесса подробно описана в работе [9].
Алгоритм продвижения дислокации в случайном массиве препятствий включает в себя комбинацию в определенной последовательности пяти элементарных актов:
1) прогибания дислокационного сегмента с нулевого напряжения до заданного;
2) расширения участка дислокации между двумя точками закрепления с начального напряжения т0 до приложенного т1;
3) проверки дислокационного узла на стабильность и (в случае невыполнения хотя бы одного из условий стабильности) прорыва элементом дислокации препятствия, определения последующей конфигурации;
4) нахождения места термической активации (стопора, который будет преодолен термически активируемым способом);
5) определения места термической активации.
При этом основу алгоритма составляет формализация первых трех пунктов.
П р о г и б а н и е д и с л о к а ц и о н н о г о с е г м е н т а. При увеличении напряжения до т1 дислокационному отрезку, закрепленному в точках А и В, представляется возможность прогибаться между этими точками в дугу окружности,радиус которой определяется напряжением т1. Процесс прогибания приводит к одному из исходов: а) прогибающийся отрезок встречает новый стопор массива - стопор С; б) отрезок беспрепятственно прогибается до равновесного радиуса кривизны; в) радиус дуги равновесного радиуса кривизны меньше половины длины отрезка АВ. Точка С находится из элементарных геометрических соотношений, основанных на теореме синусов: все вписанные в окружность углы, стороны которых проходят через две заданные точки окружности А и В, а вершина лежит по одну сторону от прямой, соединяющей эти точки, равны между собой. Поэтому стопор С выбирается так, чтобы угол АСВ был максимальным среди всех стопоров из области поиска, которая, в силу связи между напряжением и кривизной дислокационного сегмента, однозначно описывается максимально и минимально возможными углами АСВ. В случае "в" стопор А перерезается по механизму Орована.
Таким образом, на основе синтеза моделей источника Франка-Рида и скольжения бесконечных прямолинейных дислокаций через случайную сетку точечных препятствий построена модель дислокационного источника в плоскости скольжения, содержащей случайно расположенные точечные препятствия с заданным спектром прочностей. Модель позволяет в единых условиях анализировать эмиссию дислокационной петли источником и последующую эволюцию замкнутой планарной дислокационной петли в плоскости кристаллографического скольжения с барьерами для дислокаций.
Предложенная модель алгоритмизирована при активном участии к.ф.-м.н. А.В.Ушакова, за что авторы ему искренне благодарны, разработан комплекс соответствующих программ на алгоритмическом языке Фортран-77. Входными параметрами программ являются координаты источника или точки закрепления дислокационной петли, температура, относительные концентрации и прочности препятствий, характеристики кристалла, приложенное напряжение. Выходные характеристики зависят от целей исследований. Но базовую информацию представляют координаты и типы препятствий на дислокации в любой требуемой конфигурации, времена ожидания термических активаций, число пройденных стопоров при переходе из одной конфигурации в другую, число атермически пройденных стопоров после одной активации.
Литература
1. Зайцев С.И. Моделирование движения дислокаций через точечные препятствия // Дефекты в кристаллах и их моделирование на ЭВМ.- Л.: Наука, 1980.- С. 178-191.
2. Предводителев А.А. Возможности моделирования процессов, связанных с движением и размножением дислокаций в кристаллах // Динамика дислокаций. -Киев : Наукова думка, 1975.- С. 178-190.
3. Колупаева С.Н., Старенченко В.А., Попов Л.Е. Неустойчивости пластической деформации кристаллов.- Томск: Изд-во Том. ун-та,1994.- 300 с.
4. Слободской М.И., Кобытев В.С., Попов Л.Е. Зависимость средней площади, заметаемой дислокационной петлей после одной термической активации // Изв. вузов. Физика.- 1985. N3.- С. 119-120.
5. Выдашенко В.Н., Ландау А.И. Статистические характеристики конфигурации дислокации, движущейся при низких температурах // Физ. низких темпер.- 1979. Т.5, N7. - С. 794-805.
6. Форсайт Дж., Малькольм М., Моулер К. Машинные методы математических вычислений.- М.: Мир, 1980.-280 с.
7. Кнут Д. Искусство программирования для ЭВМ. Т. 2: Получисленные алгоритмы.-М.: Мир, 1977.- 924с.
8. Леман Э. Проверка статистических гипотез.- М.: Наука, 1964. - 496с.
9. Кобытев В.С., Слободской М.И., Руссиян А.А. Моделирование на ЭВМ процессов взаимодействия и скольжения дислокаций. - Томск: Изд-во Том. ун-та, 1992. - 180 с.
10. Hanson K., Morris J.W., Jr., Altintas S. Computer simulation of dislocation glide through field of point obstacles // Nucl. Met. 1976, v. 20. - P. 917-928.
11. Голосова Т.Н. Имитационное моделирование элементарных процессов размножения дислокаций в поле дискретных стопоров: Автореф. дис. ... канд. физ. - мат. наук. Томск, 1990. - 23 с.
12. Голосова Т.Н., Слободской М.И., Попов Л.Е. Моделирование источника дислокаций в поле активируемых и неактивируемых дискретных препятствий //Изв. вузов. Физика.- 1992. N10.- С.20-24.
13. Куринная Р.И.,Ганзя Л.В.,Попов Л.Е. Сопротивление расширению дислокационной петли в ГЦК металлах // Изв. вузов. Физика.- 1982. N8.- С.35-38.
14. Арнольд В.И. Математические методы классической механики.- М.: Наука,1989.-472с.