УДК 538.911, 539.32
И.Ю. Зубко, П.В. Трусов
Пермский государственный технический университет
ОПРЕДЕЛЕНИЕ УПРУГИХ ПОСТОЯННЫХ ГЦК-МОНОКРИСТАЛЛОВ С ПОМОЩЬЮ ПОТЕНЦИАЛА
^ *
МЕЖАТОМНОГО ВЗАИМОДЕИСТВИЯ
Исследуется вопрос об установлении связи между параметрами потенциала Леннарда-Джонса и упругими модулями металлов с идеальной ГЦК-решеткой. Для этого монокристалличе-ский куб подвергается различным видам деформации, находятся компоненты тензора напряжений Коши ст без априорного предположения о его симметрии. Получены аналитические выражения для коэффициентов разложений компонент тензора напряжений в степенные ряды, по ним определяются известные макроскопические линейные упругие модули. Показано, что тензор линейно-упругих свойств для металлов с ГЦК-решеткой является анизотропным и симметричным. Получена связь периода решетки ГЦК монокристалла с равновесным расстоянием для изолированной пары атомов. Получено представление второго параметра потенциала Леннарда-Джонса через период решетки и макроскопический модуль сдвига О, проведена верификация полученных соотношений. Полученные результаты могут использоваться для уточненного расчета механических свойств наночастиц.
Ключевые слова: ГЦК-монокристалла, прямое вычисление упругих модулей, нелинейная несимметричная упругость, кубическая анизотропия, параметры потенциала Леннарда-Джонса.
Введение
В связи с интенсивным развитием вычислительной техники расширяются возможности дискретных подходов описания механического поведения твердых тел, в том числе их неупругого деформирования и разрушения [1-2]. Одним из таких подходов является метод молекулярной динамики/статики, позволяющий с помощью задания потенциала межатомного взаимодействия описывать различные физикомеханические эффекты, сопровождающие деформирование твердых тел. В частности, при анализе упругого поведения твердого тела дискретный подход позволяет провести численный эксперимент по установлению особенностей упругого отклика. Например, оценить его линейность/нелинейность по отношению к накладываемым деформаци-
*
Работа выполнена при финансовой поддержке РФФИ (гранты №10-08-00156-а, № 11-01-96033-р-Урал-а).
ям, исследовать симметрийные, анизотропные свойства упругого закона, определить зависимость упругих модулей от размеров тела (в том числе - для тел микро- и нано- масштабов). Полученные зависимости могут использоваться при построении многоуровневых упругопластических моделей деформирования твердых тел, например металлических поликристаллов, мезоуровень которых отождествляется с отдельным зерном или фрагментом, в микро- и наномеханике, для уточнения моделей, используемых в механике микронеодноородных (композиционных) материалов.
Межчастичные потенциалы в дискретном подходе являются приближенным способом описания силового взаимодействия частиц материала, качественно отражающим основные свойства атомов отталкиваться на малых расстояниях и притягиваться на больших расстояниях. Получаемые количественные результаты зависят от числовых значений параметров потенциалов. В настоящее время не существует единых методик определения этих параметров для произвольных материалов. В представляемой работе предложен подход к идентификации параметров потенциала межатомного взаимодействия для металлических монокристаллов с ГЦК-структурой. Возможности подхода продемонстрированы на примере потенциала Леннарда-Джонса. С помощью этого подхода исследованы упругие свойства ГЦК-монокристаллов произвольного металла с соответствующим строением.
Для анализа упругого отклика монокристаллическое бездефектное тело в форме куба (рис. 1, а) подвергается различным видам деформации, например простого сдвига (рис. 1, б). На его гранях в деформированной конфигурации в предположении наличия только центрального характера межатомного взаимодействия определяются силы, действующие на атомы из выбранных внешних граней со стороны всех остальных атомов тела. Полученные силы делятся на площади деформированных граней, и по ним находятся компоненты тензора напряжений Коши без априорного предположения о его симметрии. Затем полученные выражения для компонент тензора напряжений раскладываются в степенные ряды по параметру, характеризующему степень деформации, например, величине сдвига у. Не зависящие от величины деформации коэффициенты при первых степенях получаемых рядов являются искомыми упругими постоянными монокристалла.
Если исследуемый материал является линейно-упругим по Коши с несимметричным тензором напряжений, удовлетворяющим обобщенному закону Гука с несимметричными мерами: о = С : Уи , о Ф от, Сум Ф с]ш, ст ф Су1к, Сук1 = Сщ, где V - оператор градиента из от-
счетной конфигурации, то полученные коэффициенты при первых степенях параметров деформации - компонент тензора дисторсии Уи -будут компонентами тензора линейно-упругих свойств С. Если материал является (нелинейно) упругим по Грину, то из разложения упругого потенциала в степенной ряд по дисторсии получается упругий закон также в виде ряда, первый член которого совпадает с правой частью записанного несимметричного обобщенного закона Г ука.
В работах [1-2] отмечается необходимость проверки условия положительной определенности тензора линейно-упругих свойств материала, компоненты которого вычисляются с помощью дискретного подхода, в частности, упоминается, что для потенциала Леннарда-Джонса это условие не выполняется. В представляемой работе проверке положительной определенности тензора С уделяется особое внимание.
Материал с кубической решеткой имеет три взаимно ортогональные оси симметрии четвертого порядка [3], совпадающие с кристаллографическими осями (см. рис. 1, а). С учетом всех симметрийных и анизотропных свойств тензор линейно-упругих свойств С имеет четыре независимые компоненты, в качестве которых выбираются компоненты с наименьшими значениями индексов в указанных кристаллографических осях: С1111, С1122, С1212, С1221. Остальные компоненты либо равны нулю, либо выражаются через них: С1133 = С3311 = С2233 = С3322 = =С =С , С =С =С =С =С =С , С =С =
С2211 С1122 , С1313 С2121 С2323 С3131 С3232 С1212 , С1331 С3113 = С2332 = С3223 = С2112 = С1221 , С2222 = = С3333 = С1111 . Ограничения на значения упругих модулей, следующие из условия положительной определенности тензора линейно-упругих свойств С ( е: С : е > 0, Уе Ф 0; е: С : е = 0 е = 0 ), получаются из неравенства
е : С : е = (еп + е22 + ез3)С1111 + 2(е11 е22 + е11 е33 + е22 е33)С1122 + (1)
+ (е12 + е21 + е13 + е31 + е23 + е32)С1212 + 2(е12 е21 + е13 е31 + е23 е32)С1221 > 0.
Из анализа этого выражения при различных еу] или с помощью критерия Сильвестра, примененного к матрице 9 х 9 с компонентами
С,у}{к/}, где каждая пара индексов {,у}, {к, 1} пробегает значения из набора пар вида {{1,1}, {1,2}, {1,3}, {2,1}, {2,2}, {2,3}, {3,1}, {3,2}, {3,3}}, выводятся следующие ограничения на компоненты:
Четвертое (двойное) неравенство (*) является нестрогим, хотя и обеспечивает строгую положительную определенность тензора С. Это связано с тем, что последние два слагаемых в неравенстве (2) в случае равенств модулей вида С1221 = ±С1212 после приведения подобных дают
строго положительную при ненулевых компонентах несимметричного тензора е сумму полных квадратов сумм или разностей:
«Несимметричный» обобщенный закон Гука в компонентах при-
В опыте на одноосное растяжение (например, при а11 Ф 0 и остальных о у = 0) может быть найден коэффициент Пуассона
нимает вид
°12 — С1221(^и)12 ^ С1212(^и)21’ 013 = С1221(^и)13 + С1212(^и)31’ °21 — С1212(^и)12 ^ С1221(^и)21’
(3)
°23 — С1221(^и)31 + С1212(^и)32’ °31 — С1212(^и)13 ^ С1221(^и)31’ °32 — С1212(^и)23 + С1221(^и)32’
033 — СЦ22(^и)11 + С1122(^и)22 + С1111(^и)33'
Это выражение монотонно возрастает с ростом C1122 / С1Ш в интервале значений, обусловленном неравенством (2). С использованием условий (2) получается ожидаемый результат в виде ограничений на значения коэффициента Пуассона:
-1 < v < 1/2 . (5)
1. Простой сдвиг монокристалла
Пусть n - единичная нормаль к плоскости сдвига, m - единичный вектор, задающий направление сдвига. Закон движения в лагранжевой форме имеет вид: r = R + у(n• R)m, где R - радиус-вектор центров
атомов в отсчетной конфигурации, а r - радиус-вектор атомов в текущей конфигурации. Тогда деформационный градиент F = (Vr)T , описывающий сдвиг на величину у, имеет вид F = E + ymn . Рассмотрим верхнюю и правую боковую грани. Площадь верхней грани при описанном сдвиге не меняется и равна N2а2, так как на каждый атом на грани приходится ее часть площадью а2, где а - параметр решетки, N - количество атомов вдоль ребра куба (на рис. 1 N = 4). Площадь деформируемой грани меняется согласно соотношению d s = J(N • U 2 • N)12 d S , где N - единичный вектор нормали к недефор-мированной грани, U 2 = F 1 • F-T , F 1 = E-ymn, J = detF, dS - площадь грани до деформирования, d s - площадь грани после деформирования [4]. Для правой боковой грани получим площадь после сдвига,
равную N2а2ф + у2 .
Выражение для суммарной силы, действующей на любую из граней, для произвольной силы f (|r|) r / |r| взаимодействия пары
атомов, соединенных вектором r, могут быть найдены аналитически с помощью любого пакета символьных вычислений, например пакета Wolfram Research «Mathematical. Для произвольных значений межатомного расстояния а может оказаться, что в состоянии до приложения сдвига суммарные силы на гранях исходного куба отличны от ноля, т.е. в качестве начальной конфигурации берется не естественное состояние (ненапряженное и недеформированное), а некоторое равноосно деформированное состояние. Для получения
естественного состояния необходимо решать задачу о выборе естественного периода решетки. Эта проблема учитывается далее при нахождении упругих модулей.
Рис. 1. Монокристаллическое тело из 27 элементарных ячеек ГЦК-решетки: а - до сдвига, б - после простого сдвига в плоскости X1OX3 в направлении X1. Вдоль осей координат - межатомные расстояния
Нормаль n к деформированной грани определяется выражением n = (N• U-2 • N)-1/2F -T • N [4]. Например, для правой боковой грани с нормалью NR = {1,0,0} в недеформированной конфигурации после рассматриваемого простого сдвига получим нормаль в текущей конфигурации nR = Ц/^Д + у2, -у / -^1 + у2 ,о|. Для нахождения матрицы компонент тензора напряжений Коши а у в приведенном на рис. 1 базисе кристаллографической системы координат используется соотношение Коши n • о|г = tn |г . Это соотношение, записанное для правой боковой, верхней и передней граней, в рассматриваемом случае позволяет найти все 9 компонент тензора напряжений Коши: ап(у) = у( tUp ) +
+V1 + У2 (tR )l, а12(У) = y(^p )2 +V1 + у2 (tR ^2 , а21(У) = (Ч )1, а22(У) =
= (tUp )2, а13 = (tUp )3 + V1 + у (tR )з , а23 = (tUp )3 , а31 =(^Fr )1, а32 =(^Fr )2 ,
а33 =( tFr )3, где индекс «R» (Right) относится к правой боковой грани,
индекс «Up» (Upper) - к верхней, индекс «Fr» (Front) - к передней грани. В частности, полученные выражения дают возможность сравнивать недиагональные компоненты тензора напряжений а12(у) и а21(у). Для
определения упругих модулей функции а12(у) и а21(у) раскладываются в степенные ряды в окрестности ноля ( у = 0), в которых коэффициенты при первых степенях будем называть модулями сдвига С12 и С21: а12(у) = а12(0) + С12у + 2а"2(у)|у_0 У2 + ••• Значения компонент тензора
напряжений Коши в начальной конфигурации а12(0) и а21(0) не оказывают влияния на значения исследуемых модулей. Для линейноупругого материала найденные таким образом касательные (упругие)
« П> а12(у) — СТ12(0)
модули будут совпадать с секущими модулями ц2 - —12------------------—
дули. Для аналитического определения введенных коэффициентов использовался пакет «МаШешайса 7.0», в котором заложены довольно эффективные алгоритмы символьных вычислений.
Закон Гука (1) для рассматриваемого простого сдвига принимает вид
То есть введенные модули сдвига С12 и С21 связаны с компонентами тензора С как
Выберем в качестве потенциала взаимодействия атомов (двух-
|г| _ г , где г - вектор, соединяющий центры двух атомов. Тогда модуль силы взаимодействия / (г) _ —ф'(г) . Используя конкретный вид силы,
можно в явном виде получить набор коэффициентов разложения в ряды Тейлора компонент а12(у) и а21(у) тензора напряжений Коши в окрестности ноля (у_ 0). Заметим, что при аргументе г _а сила взаимодействия атомов обращается в ноль, т.е. параметр а есть равновесное расстояние для изолированной пары атомов, но а Ф а, так как
у
„ Г* - а21(у) — а21(0)
и ^21 -----------------------
у
. Далее рассматриваются только касательные мо-
а12 ^1212у, а21 _ ^1221у,
(6)
(7)
частичный) потенциал Леннарда-Джонса
параметр решетки а определяется из условия равновесия большого количества атомов. Действительно, даже для одномерной цепочки из 3 атомов крайние атомы приблизятся к центральному атому, и равновесное расстояние а в такой системе будет отличаться от величины а.
Для ГЦК-решетки, начиная с количества атомов N _ 2 , были получены коэффициенты разложения в ряд Тейлора при первых девяти степенях (часть результатов приведена в табл. 1-4). Для коэффициентов степенного ряда компоненты а12 вводятся обозначения:
а12(у) _ а12(0) + в12у + у2 + О™ у3 +... Все коэффициенты имеют одина-
6 ^ ” 6 ^ ••
ковое строение вк _ а (С1 • а + С^ а6 )р / а15 и отличаются числовыми
15
множителями С1 • (N) и С2 • (N), зависящими от числа атомов N.
Таблица 1
Коэффициенты при 1-й степени у ГЦК-решетки
N Коэффициент в12 Коэффициент в21
5 а6 (—570,51а6 +11061,12а6)р /а15 а6 (—570,51а6 +11061,12а6)р/а15
9 а6 (—673,15а6 +12901,44а6 )р/ а15 а6 (—673,15а6 +12901,44а6)р/а15
12 а6 (—706,64а6 +13501,53а6)р/ а15 а6 (—706,64а6 +13501,53 а6)р/ а15
Таблица 2
Коэффициенты при 3-й степени у ГЦК-решетки
N Коэффициент в ш12 Коэффициент в Ш21
5 а6 (—1132,52 а6 + 77545,23 а6 )р / а15 а6 (—1171,63а6 + 77571,52а6)р/а15
9 а6 (—1302,19а6 + 90403,51а6)р/а15 а6 (—1361,01а6 + 90441,51а6)р/а15
12 а6 (—1357,20а6 + 94596,15а6)р/а15 а6 (—1422,81а6 + 94638,12а6)р/а15
Таблица 3
Коэффициенты при 5-й степени у ГЦК-решетки
N Коэффициент 0У12 Коэффициент вУ21
5 а6 (2649,56а6 +14919,42а6)р /а15 а6 (2634,43 а6 +14993,94а6 )р / а15
9 а6 (3074,33 а6 +17424,02 а6 )р / а15 а6 (3052,91а6 +17531,50а6)р/а15
12 а6 (3212,95 а6 +18239,82 а6 )р / а15 а6 (3189,45а6 + 18358,71а6)р/а15
Оказалось, что коэффициенты при первых степенях у при разложении компонент а12(у) и а21(у) совпадают. Коэффициенты при четных степенях на 15-18 порядков меньше коэффициентов при первых степенях и ими можно пренебречь. Коэффициенты при третьих, пятых и более высоких нечетных степенях у для разложений в ряды Тейлора а12(у) и а21(у) оказались различными. Но при малых значениях сдвига
(до величины у = 0,1) члены разложения с третьей и выше степенями у дают пренебрежимо малый вклад в выражения для компонент а12 (у) и а21 (у) . Дадим оценку вклада этих членов, принимая для простоты, что а «а (см. табл. 4).
Таблица 4
Коэффициенты при нечетных степенях у ГЦК-решетки
N Оатаъ/ Р аъОш12 / Р Р/ Р а30У12/ Р Р/ > ^ Р/ > Р/ >
5 10490,6 76412,7 76399,9 17569,0 17628,3 -150942,0 -150961,0
9 12228,3 89101,3 89080,5 20498,3 20584,4 -176073,0 -176101,0
12 12794,7 93238,9 93215,3 21452,8 21548,1 -184266,0 -184297,0
Для сдвига на величину у = 0,1 и числа атомов N = 12 на ребре исходного куба от слагаемого третьего порядка получим добавку 932,389Р / а3 к значению компоненты а12 и добавку 932,153Р / а3 к значению компоненты а21. Относительное различие в полученных значениях составит примерно 1,72 х10-3%. Для величины сдвига у = 0,05 различие будет значительно меньшим. Таким образом, отклонение тензора напряжений от симметричного вида для монокристалла с ГЦК-решеткой в условиях малых упругих деформаций мало, определяется членами ряда с более высокими степенями параметра сдвига у и зависит от величины у. Упругий закон при этом оказывается нелинейным.
Для числовых коэффициентов из табл. 1 найдены аппроксимирующие функции и их предельные значения для числа атомов N ^ да . Например, графическое представление зависимостей от числа атомов N первого и второго числовых коэффициентов Сщ (^ и Сщг (^,
стоящих в выражении для модуля С12 из табл. 1
- а6 (С,Ч:а‘ + ещ,а‘)р/ а15, приведены на рис. 2. Для С1С.ДN и
С2сь (N было проверено, что они не сводятся к таким, не имеющим горизонтальной асимптоты функциям, как логарифмическая и степенная функция с меньшим единицы положительным показателем. Окончательно зависимости искались в классе функций /(х) = а(х - х0 )к + Ь . Методом наименьших квадратов получены аппроксимирующие выражения С^(N) = 1343,13(N + 0,47)-101 -811,86 (при N ^да коэффициент С1б2
стремится к значению -811,86) и С2ва(Щ = -23608,90(N + 0,41)-1,008 + +15367,82 (при N ^да этот коэффициент стремится к значению 15367,82). Показатели степени близки к -1, то есть кривые на рис.2 являются ветвями гипербол с центрами в точке с абсциссой, примерно равной -0,45. Итак, упругий модуль Б12 представляется в виде
+
(8)
а6 ((1343,13( N + 0,47)-ш - 811,86) а + (-23608,90( N + 0,41)-1’008 +15367,82) а6 )р / а15, а его предельное при N ^ да (макроскопическое) значение находится как С12да = а6 (-811,86а6 + 15367,82а6 )р / а15. (9)
С1(С12)
■ N
Рис. 2. Зависимость числовых коэффициентов С10 (N и С20 (N) от числа атомов N
Для модулей более высокого порядка получены аналогичные выражения, причем все числовые коэффициенты, вычисленные для различных значений N ложатся на ветви подобных гипербол с центрами в точках с абсциссами из интервала (-0,47; -0,40). Выражения для остальных модулей в статье не приведены.
В упрощенном случае, когда для получения оценок принимается приближенное равенство периода решетки и параметра потенциала а «а , графики полученных зависимостей недиагональных компонент тензора напряжений а12(у) и а21(у) в виде рядов с предельными
(N ^ да) значениями числовых множителей приведены на рис. 3, а. Относительное несовпадение этих компонент в зависимости от величины сдвига показано на рис. 3, б. Как видно из полученных результатов, относительное различие между недиагональными компонентами тензора напряжений а12(у) и а21(у) в опыте на простой сдвиг ГЦК-
монокристалла при упругом деформировании составляет менее одной сотой процента и им можно пренебречь. Поэтому можно говорить только об одном модуле сдвига ГЦК-материалов Б = Б12 = Б21. Другие
анизотропные модули СШ1 и СШ2 определяются в опыте на чистое
растяжение-сжатие.
0,008
/ \ 0,006 ■ \
\ 0,004
(1002 ]
-0,4
-0,2
0,2
0,4
Рис. 3. Зависимость компонент ст12(у) и ст21(у) от величины сдвига у для ГЦК-решетки: а - нормированные значения компонент, б - относительное несовпадение значений компонент ((ст12 -ст21)/ ст12 )х100%
Сами зависимости компонент а12(у) и а21(у) от величины сдвига у для ГЦК-решетки являются нелинейными, а полученные разложения в степенные ряды справедливы только в окрестности ноля. Т очное выражение для этих компонент, применимое для произвольных значений параметра сдвига у, хотя и выведено, но не может быть представлено в полном виде. Использование пакета символьных вычислений позволяет представить зависимость для него в графическом виде (например, на рис. 4 приведен график для случая N = 12).
Отметим еще одну особенность проводимых вычислений. При определении напряжений на гранях исходного недеформированного куба все касательные напряжения оказываются нулевыми, а нормальные совпадают между собой и выражаются через параметры потенциала а, р и параметр решетки а. Можно потребовать, чтобы эти напряжения также были равны нулю, что будет соответствовать получению естественного ненапряженного начального состояния. Из этого требования выводится однозначная связь параметров решетки а и потенциала а (рис. 5). То есть в пределе при N ^ да с использованием потенциала межатомного взаимодействия Леннарда-Джонса все упругие модули будут выражены через два параметра а и р.
Рис. 4. Графики функций ст12(у) Рис. 5. Зависимость параметра решетки а
и ст21 (у) (неразличимы), N = 12 от параметра потенциала а
Леннарда-Джонса
Приведенная на рис. 5 зависимость описывается соотношением
а = 0,0104 (N + 0,0834 V1’03 +1,3918, (10)
а
т.е. для ГЦК-решетки и потенциала Леннарда-Джонса в пределе ( N ^ да ) получим
а, = 1,3918 а
(11)
То, что период ГЦК-решетки а, оказался больше, чем параметр потенциала а равновесного расстояния для изолированной пары атомов, объясняется строением ГЦК-кристалла, на гранях элементарных ячеек которого расположены дополнительные атомы, «раздвигающие» атомы, лежащие на ребрах ячеек. Для простой кубической решетки период а, будет меньше, чем параметр а - в аналогичных расчетах получено значение а, = 0,9561а , для ОЦК-решетки получена связь а, = 1,1251а .
Исследуем вопрос об однородности поля напряжений в изучаемом объеме атомов. Рассматривается начальная конфигурация моно-кристаллического куба для случаев с различным числом атомов N на ребре. Период решетки а в зависимости от N берется различным согласно соотношению (10), чтобы обеспечить нулевые напряжения на поверхности куба. Для полученной (отсчетной) конфигурации строится набор сечений по атомным плоскостям, параллельным граням куба. Часть атомов, отделенных от тела проведенным сечением отбрасывается, а ее влияние заменяется эквивалентными ей компенсирующими силами. В каждом таком сечении вычисляется вектор напряжений. Для этого находится сумма всех сил, действовавших на атомы сечения со стороны всех отброшенных атомов тела, и делится на площадь сечения. Показано, что получаемый вектор напряжений в рассматриваемой (отсчетной) конфигурации всегда направлен строго по нормали к сечению. Получена зависимость единственной ненулевой компоненты вектора напряжений от положения сечения в объеме куба. Для примера приводятся результаты для правой боковой грани куба и вычисляется зависимость а11 от к, к = 1, N. Оказалось, что при постоянном периоде
решетки однородность (отличие от ноля) напряжений нарушается только в тонком поверхностном слое, толщина которого составляет порядка 2-3 межатомных расстояний (рис. 6, а), и в нем возникают сжимающие напряжения. Это говорит о том, что в этом поверхностном слое атомные плоскости стремятся сблизиться и образовать пленку с большей плотностью, чем у материала в объеме. Максимальное значение найденной компоненты тензора напряжений, достигаемое на второй от поверхности атомной плоскости, зависит от числа атомов N (рис. 6, б) и описывается зависимостью
Ч N) -
-0,463 + -
V
(-0,179 + N)
1,001
а
Задача отыскания неоднородного вблизи свободной поверхности распределения атомных слоев и определения энергии поверхностного натяжения не решалась, эти вопросы будут рассмотрены в последующих публикациях. Заметим, что поверхностная энергия будет зависеть от размеров тела, что косвенно подтверждается приведенной зависимо-
_ max / л т\
стью ай (N).
Рис. 6. Изменение напряжений по длине образца: а - распределение значений ст11 по длине образца при N = 50; б - зависимость пикового значения аи = ст11,
от числа атомов N
При реализации простого сдвига распределение значений нормальных компонент ай тензора напряжений Коши не меняется. На рис. 7 приведена зависимость значений касательных напряжений а12, определяемых на внешней грани тела в опыте на простой сдвиг на величину у = 105 от числа атомных слоев к, взаимодействие с которыми учитывается при расчете а12 (к откладывается от поверхности в объем тела). Таким образом, получаемое напряженное состояние можно считать однородным.
Итак, параметр потенциала а Леннарда-Джонса для моноатом-ных кристаллов с бездефектной ГЦК-решеткой выражается как
а = а/1,3918. (12)
Модуль сдвига при этом О = 66,4548 Р / а3 согласно (9), следовательно, второй параметр потенциала Леннарда-Джонса р находится по формуле
р=
66,4548
= 0,0150 в а3 = 0,05581 в а,
(13)
Рис. 7. Изменение касательных напряжений по длине образца при у = 10-5, N = 50 : а - зависимость ст12 от числа слоев к, используемых в расчете, б - относительное отклонение ст12 (к) от значения ст12 (N -1)
Для меди ( а„Си = 0,3615 х109 (м), ССи = 4,55 х1010 (Па)) с использованием формул (12) и (13) получим значения параметров потенциала Леннарда-Джонса:
аСи = 2,5974 х10-10 (м), рСи = 1,1996 х10-20 (Н-м). (14)
Приведенные соотношения позволяют найти зависимость упругих модулей кристалла от числа атомов N на ребре исследуемого объема или при известном периоде решетки от размеров образца Ь = N а„. Графическое представление зависимости модуля сдвига и периода решетки от числа атомов N для меди приведено на рис. 8.
Рис. 8. Зависимость от числа атомов N на ребре куба: а - модуля сдвига меди,
б - периода решетки меди
Вычисления с найденными параметрами потенциала показывают, что модуль сдвига для образца меди с размером 1 см будет составлять = 4,54х1010 (Па), для образца с размером 1 мкм -
= 4,53х1010 (Па), с размером 100 нм - = 4,51х1010 (Па),
с размером 10 нм - = 4,28 х1010 (Па). Оказалось, что с уменьше-
нием размеров монокристалла его модуль сдвига также уменьшается, и существует такой «переходный» размер образца (рис. 8, а), что при дальнейшем сокращении размеров тела падение модуля происходит очень быстро. С уменьшением размеров тела равновесное межатомное расстояние растет, следовательно, уменьшается плотность частиц
с ГЦК-решеткой ( р(N = 4та / (а(Ы))3, где та - масса атома).
2. Чистое растяжение-сжатие монокристалла
Далее исследуется поведение материала в опыте на чистое растяжение-сжатие исходного объема с ГЦК-решеткой (рис. 9) при сохранении размеров в двух других направлениях (ортогональных оси растяжения-сжатия). Пусть т - единичный вектор, задающий направление оси растяжения-сжатия (для деформации, показанной на рис. 8, т = {1,0,0}). Деформационный градиент, описывающий соответствующую деформацию с кратностью удлинения X, имеет вид Г = Е + (X- 1)тт.
Рассмотрим верхнюю, переднюю и правую боковую грани куба, задаваемые соответственно нормалями ^р = {0,1,0}, = {0,0,1},
N к = {1,0,0} в отсчетной конфигурации. Площадь правой боковой грани при описанном деформировании не меняется и равна Nа2. Площади деформируемых (верхней и передней) граней равны между собой и меняются согласно соотношению = 3^ • и 2 • ^1/2ё£ [4], где N -единичный вектор нормали к грани в отсчетной конфигурации, и 2 = Г-1 • Г-т, Г-1 = Е + (X-1 - 1)тт , 3 = ёе! Г = Х , - площадь грани до деформирования, - площадь грани после деформирования.
Для верхней и передней граней получим одинаковые значения их площади ё£ир = XN2а2. Нормали к граням при растяжении-сжатии не меняются: п„р = Nир = {0,10}, пРг = NРг = {0, а ^ пк = Nк = {1,0,0} .
а б
Рис. 9. Монокристаллическое тело из 27 элементарных ячеек ГЦК-решетки: а - после растяжения вдоль ОХ^; б - после сжатия вдоль ОХь Вдоль осей координат - межатомные расстояния
Для нахождения матрицы компонент тензора напряжений Коши о у в приведенном на рис. 9 базисе кристаллографической системы координат используется соотношение Коши п • о|г = 1Й |г . Записывая это
соотношение для верхней, передней и правой боковой граней, в рассматриваемом случае получим выражения для компонент этого тензора: ап(у) = (^ )р Му) = (^ ^2, °13 =(¿к ^ о21(у) = (^ )1,
°22(у) = (¿Ир )2 , °23 =(¿Ир )3 , °31(у) = (^г )1 , °32(у) = (^г )2 , °33 = (^г )3 .
Для определения упругих модулей функции ап(Х), о22(Х) и о33(Х)
раскладываются в степенные ряды по параметру X в окрестности точки Х = 1. Коэффициент при первой степени в разложении о11(Х) будем
обозначать как модуль растяжения-сжатия Е11, а коэффициенты при первых степенях в о22 (X) и о33 (X) разложении - как модули растяжения-сжатия Е22 и Е33. При этом значения компонент тензора напряжений в начальной конфигурации не будут влиять на значения введенных модулей. Для линейно-упругого материала найденные таким образом касательные (упругие) модули будут совпадать с секущими модулями. Из компонентной формы записи обобщенного несимметричного закона Гука (3) для рассматриваемого опыта на чистое растяжение-сжатие получим
о11 = С1111(Х-1),
< °22 = °33 = ^1122 (Х— 1),
°12 = °21 = °13 = °31 = °23 = °32 = 0
Следовательно, введенные модули Е11, Е22 и Е33 связаны с компонентами тензора линейно-упругих свойств С как
Е11 = С1111 , Е22 = Е33 = С1122 . (16)
Для ГЦК-решетки при использовании потенциала Леннарда-Джонса в опыте на чистое растяжение-сжатие были получены коэффициенты разложения в ряд Тейлора при нескольких степенях (часть результатов приведена в табл. 5-9). Для коэффициентов степенного ряда (упругих модулей) рассматриваемых компонент тензора напряжений вводятся обозначения вида о11(Х) = о11(1) + Е11(Х — 1) + Е“(Х — 1)2 + Е11|1(Х — 1)3 +..., где
числа Ек, К = II, III,... являются упругими модулями К-го порядка.
Таблица 5
Коэффициенты при 1 -й степени X -1
N Коэффициент Е11 Коэффициент Е22 = Е33
5 а6 (—1894,99 а6 + 24470,32 а6 )р / а15 а6 (—1784,91а6 + 19962,91а6 )р/а15
9 а6 (—2206,63а6 + 31275,81а6)р/а15 а6 (—2094,32 а6 + 23279,43 а6 )р / а15
12 а6 (—2308,48 а6 + 32723,31а6 )р / а15 а6 (—2195,65а6 + 24613,51а6)р/а15
Таблица 6
Коэффициенты при 2-й степени X — 1
N Коэффициент Е" Коэффициент Е22 = Е3П3
5 а6 (2977,55 а6 — 79354,93 а6 )р / а15 а6 (4125,23 а6 — 74242,62 а6 )р / а15
9 а6 (3442,55 а6 — 90306,31 а6 )р / а15 а6 (4845,15 а6 — 86581,52 а6 )р / а15
12 а6 (3594,51а6 — 96658,21а6)р/а15 а6 (5080,15 а6 — 90605,32 а6 )р / а15
Коэффициенты при 3-й степени X-1
N Коэффициент Е™ Коэффициент е£ = Е™
5 а6 (-3148,25 а6 +148743,01 а6 )Р / а15 а6 (-7112,22а6 +198365,04а6)р/а15
9 а6 (-3564,82а6 +172791,02а6)Р /а15 а6 (-8369,09а6 + 231345,02 а6 )р / а15
12 а6 (-3702,02а6 +180639,06а6)Р/а15 а6 (-8781,11 а6 + 242101,01 а6 )р / а15
Таблица 8
Коэффициенты разложения для Е^ ГЦК-решетки
N Е11а3 / р а3 Е“/ Р са а3Е" / Р а3 Е^/ Р а3 Е"/ Р а3 Е"1/ Р
5 24942,4 -76377,4 145594,0 -188333,0 179380,0 -213911,0 611768,1
9 29069,2 -88960,1 169226,1 -217272,3 200850,4 -223702,1 642140,4
12 30414,8 -93063,7 176937,3 -226733,2 207937,1 -227155,2 652755,1
Таблица 9
Коэффициенты разложения для Е^ = Е3М ГЦК-решетки
N ^33 а3 е22/ Р са а3Е" / Р а3 Е2у2/ Р а3 Е"/ Р са > £
5 18178,1 -70117,3 191252,0 -390758,1 623608,2 -789652,0 781987,2
9 21185,0 -81736,4 222976,3 -455631,1 727249,2 -921113,0 912620,1
12 22165,2 -85524,3 233320,1 -476784,0 761046,2 -963988,1 955235,0
Полученные результаты согласуются с видом закона (11) и условием положительной определенности тензора С для любых металлов с ГЦК-решеткой. У коэффициентов при последующих членах разложения в степенной ряд продолжается чередование знаков, причем при нечетных степенях коэффициенты положительны. Были построены аппроксимирующие функции для найденных числовых множителей как функций от количества атомов N на ребре куба (приведена часть зависимостей):
Е11N = а6 ((3998,42(N + 0,42)-1,005 - 2626,65) а6 +
+ (-56494,62( N + 0,39) 1>004 + 37232,60 )а6 )р / а15,
Е11 „ - а6 (-2626,65а6 + 37232,60а6)Р/а15. (17)
Е22 N = Е33N = а6 ((3991,07(N + 0,43)-1’005 - 2512,43) а6 +
+ (-42206,04( N + 0,39) 1>004 + 27729,90) а6 )Р / а15,
Е22 „ = Е33 „ =а6 (-2512,43а.6 + 27729,90а6 )Р / а^5. (18)
Все найденные зависимости имеют практически равные показатели степени к « -1, то есть являются ветвями гипербол, центры которых имеют почти одинаковое значение абсциссы х0 » -0,4 .
Для выведенного ранее соотношения а. = 1,3918 а (11) получена зависимость «макроскопических» компонент а11, а22 от кратности удлинения X, представленная графически на рис. 10. Для произвольной величины растяжения (рис. 10, а, правая ветвь) наблюдается падение напряжения, связанное с разрушением образца вследствие такого расхождения атомных слоев, что сила взаимодействия между ними стремится к нулю. При сжатии, напротив, сила отталкивания слоев неограниченно возрастает по модулю из-за невозможности «бесконечно близко» сдвинуть атомы.
Макроскопические модули Е11 и Е22 = Е33 при найденной для
ГЦК-материалов связи а. = 1,3918 а выражаются как
Е11 = 127,36 Р / а3, Е22 = Е33 = 66,47 Р / а3. (19)
В частности, для меди при полученных выше значениях параметров потенциала а = 2,5974 х10-1° (м), Р = 1,1996 х10-2° (Н-м) следует, что
ЕС = 8,719 х 1010 (Па), Е2с2и = Е3с3и = 4,551 х 1010 (Па). (20)
В литературе имеются значения Е1с,и = 8,4 х 1010 (Па) для сверхчистой литой меди, Е^“ = 11,0 х1010 (Па) для прокатанной меди, Е1с,и = 13,0 х1010 (Па) для холоднотянутой меди. Полученное значение Е^“ из (20) довольно точно попадает в интервал приведенных экспериментальных значений макроскопического модуля Юнга меди (отличие для литой меди составляет порядка 3,5 %). Найденные значения упругих модулей (19), и в частности, для меди (20), удовлетворяют условию (2) положительной определенности тензора С.
Рис. 10. Зависимость безразмерных комплексов ст11/22а3 / Р от кратности удлинения X: а - для произвольных значений X; б - для малых значений разности X -1, сплошная линия для компоненты о11а3 / Р , пунктирная - для компоненты о22а3 / Р
Коэффициент Пуассона получился одинаковым для всех материалов с ГЦК-решеткой и согласно (4) и (19)
V — 0,343. (21)
Найденное значение коэффициента Пуассона для металлических монокристаллов с ГЦК-решеткой (21) также соответствует условию положительной определенности тензора линейно упругих свойств С в виде (5). Для меди экспериментально определенное значение коэффициента Пуассона VСи — 0,35, для алюминия - vAl — 0,34. Результаты
применения предложенной в работе методики теоретического расчета величины этого коэффициента дают совпадение с его экспериментальными значениями с точностью 2 %.
Таким образом, верификация методики дала очень высокую точность расчетных значений упругих модулей в сравнении с экспериментальными данными. Показано, что упругий закон для тел микро- и наноразмеров является анизотропным, причем свойства материала при уменьшении размеров изменяются - значения упругих модулей и плотности также уменьшаются. Наличие в реальных телах дефектов типа дислокаций, вообще говоря, должно нарушать симметрию кристалла и, следовательно, симметрию упругого отклика, т.е. вызывать отклонение тензора напряжений от симметричного вида. Однако для тел микро- и наноразмеров дислокации выходят на поверхность, а внутренняя область остается бездефектной. Для макроскопических
монокристаллов с равномерным распределением дислокаций по системам скольжения отклонение от симметрии будет нивелироваться за счет суммирования межатомных сил по большому объему. Изменить симметрийные свойства монокристалла может образовавшаяся в нем дислокационная структура, но этот вопрос лежит за рамками представленной работы.
Выводы
В работе предложена методика идентификации параметров потенциала межатомного взаимодействия Леннарда-Джонса по макроскопическому упругому модулю сдвига и макроскопическому периоду решетки чистых бездефектных металлических ГЦК-монокристаллов. Проверена однородность напряжений по объему исследуемого образца. Полученные точные выражения для компонент тензора напряжений Коши, представленные также в виде степенных рядов, показывают, что тензор линейно упругих свойств таких материалов является анизотропным и симметричным, отклонение тензора напряжений от несимметричного вида, возникающее за счет вклада нелинейных слагаемых, мало. Получена зависимость модуля сдвига от размеров образца, причем показано, что чем меньше размеры тела, тем меньше модуль сдвига. Плотность тела также уменьшается при уменьшении размеров тела за счет роста равновесного межатомного расстояния. Модуль Юнга, вычисленный согласно предложенной методике с помощью параметров потенциала Леннарда-Джонса, идентифицированных по макроскопическому модулю сдвига для меди и периоду ее кристаллической решетки, экспериментально определенный при нормальной температуре, совпадает с известным из макроскопических экспериментов значением модуля Юнга меди с точностью 3,5 %. Коэффициент Пуассона V оказался одинаковым для всех ГЦК-монокристаллов и равным 0,343, что с точностью 2 % совпадает со значением коэффициента Пуассона меди и алюминия. Все упругие модули, теоретически определенные в работе с помощью потенциала межатомного взаимодействия Леннарда-Джонса, удовлетворяют полученному для анизотропного материала с кубической симметрией условию положительной определенности тензора линейно-упругих свойств.
Библиографический список
1. Кривцов А.М. Деформирование и разрушение твердых тел с микроструктурой. - М.: Физматлит, 2007. - 304 с.
2. Кривцов А.М. Упругие свойства одноатомных и двухатомных кристаллов: учеб. пособие. - СПб.: Изд-во С.-Петерб. гос. политехн. ун-та, 2010. - 144 с.
3. Черных К.Ф. Введение в анизотропную упругость. - М.: Наука, 1988. - 190 с.
4. Поздеев А.А., Трусов П.В., Няшин Ю.И. Большие упругопластические деформации: теория, алгоритмы, приложения. - М.: Наука. 1986. - 232 с.
Получено 21.03.2011