|
Сингония
|
Характеристики элементарной ячейки
|
Принятые правила расположения осей
|
|
Триклинная
|
a ≠ b ≠ c; б ≠ в≠ г ≠ 90°
|
По рёбрам кристалла c<a <b
|
|
Моноклинная
|
a ≠ b ≠ c; б = г =90°≠ в
|
Ось y
вдоль оси 2 или вдоль нормали к плоскости m
|
|
Ромбическая
|
Оси x, y, zвдоль трёх взаимно перпендикулярных осей 2 или вдоль
нормалей к плоскостям m
|
|
Тетрагональная
|
a = b ≠ c; б = в= г = 90°
|
Ось 4 или 4¯- ось z. Остальные - в плоскости xy.
|
|
Гексагональная
|
a = b ≠ c; б = в = 90°, г=120°
|
Ось 6 или 6¯ (ось 3 или 3¯ в тригональной сингонии) ось z. Остальные в плоскости xy.
|
|
Кубическая
|
a = b = c; б = в=г = 90°
|
Оси x, y, zпараллельны трём взаимно перпендикулярным осям 4,
или 4¯ или 2.
|
Положение любой точки в кристалле обычно задается, радиус-вектором
R=ma+nb+pc (4)
где m,n,p - координаты точки в выбранной системе координат, единицами
длины в которой являются модули примитивных трансляций a, b, c. Понятно, что
координаты узлов примитивной решетки целые числа. Чтобы указать на определенный
узел его координаты заключают в двойные квадратные скобки [[mnp]].
1.2 Функции, описывающие физические величины
в кристалле.
Трансляционная симметрия кристаллов является наиболее общим свойством
кристаллов. При трансляции кристалла каждая из точек внутри кристалла занимает
место ей эквивалентное. Это значит, что во всех точках, отстоящих от данной на
любой вектор трансляции T, все физические величины во всех этих точках имеют
такие же значения, как и в исходной. Если избранная точка определена вектором
r, а значение некоторой физической величины (электронной плотности,
электрического потенциала, напряженности электрического или магнитного полей,
вектор-потециала и др.) в той точке есть f(r), то и во всех точках r+T значение
этой величины точно такое же, т.е.
f(r+T) = f(r) (1)
но T =
ma + nb+ pc,
где m,n,p - целые числа, а a,b,c - примитивные трансляции, следовательно,
функции f, описывающие физические величины в кристалле, периодические, и их
периодами являются примитивные трансляции (обозначим их для краткости ai: a1=
a, a2= b, a3 = c).
Но всякую периодическую функцию можно разложить в ряд Фурье (в нашем
случае трехмерный). Пусть ai ортогональны между собой, тогда r=ix+jy+kz,
a1=ia1, a2=ja2, a3=ka3, и
(2)
здесь
n1, n2, n3 -всевозможные целые числа. Последнее выражение можно записать в виде
.(3)
Показатель
экспоненты в (3) не что иное, как скалярное произведение двух векторов, именно
вектора r и вектора
g = in1(2p/a1)+jn2(2p/a2)+kn3(2p/a3) (4)
Целые
числа n1,n2,n3 , принимая все возможные значения, образуют бесконечное
множество векторов g
Если
вектор r определяет какой-нибудь узел решетки Браве, т.е. совпадает с вектором
трансляции T, который в данном случае имеет вид T = i a1+ j a2 + k a3, то
скалярное произведение
gT=2
(a1t1n1/a1+a2t2n2/a2+a3t3n3/a3)=2p(t1n1+t
2n2+t3n3)=2pN (5)
Здесь
N целое, т.к. ni и ti целые, и, следовательно 
,
Отсюда
следует, что
(6)
т.е.
для всех точек, отстоящих от данной на любой вектор T, значение физической
величины f одно и то же, поэтому они и эквивалентны.
Введем
обозначения
тогда
g=b1n1+b2n2+b3n3 -вектор, имеющий размерность обратной длины. Из определений
следует, что 
. Обобщая это соотношение на случай, когда ai образуют
косоугольную элементарную ячейку, можно показать, что и в этом случае
показатель степени экспоненты в разложении функций в ряд Фурье есть скалярное
произведение векторов r=a1x1+a2x2+a3x3, где xi- координаты точек в данной косоугольной
системе координат, и вектора g=b1n1+b2n2+b3n3, где
единичные
векторы в пространстве, в котором определяются векторы g.
Если
построить все эти векторы из одной точки, изображая их концы точками, то
получим систему точек, которые образуют пространственную решетку. Векторы g
являются в этой решетке трансляциями, а bi -примитивными трансляциями в ней.
Элементарная ячейка, построенная на них, имеет объём равный
(b1[b2b3])=8p3/(a1 [a2a3])
= 8p3/Vэя, (7)
здесь,
как видно, Vэя - объём элементарной ячейки в кристаллической решетке. А
решетка, образованная концами векторов g, называется обратной решеткой.
1.3 Колебания линейного многоатомного кристалла
При k , стремящемся к бесконечности, l стремится к нулю, что означает, что
приближение сплошного тела непригодно для нахождения колебаний кристалла. Пусть
имеется одномерный кристалл, составленный из одинаковых атомов с массой M, и постоянной решетки a. поло жим также, что взаимодействие
атомов подчиняется закону Гука, а именно, силы взаимодействия пропорциональны
изменению расстояния между ними (в первом приближении учтем взаимодействие
только между наиближайшими соседями). Перенумеруем атомы индексом n, а смещения атомов из положения равновесия
обозначим un. Координата атома с номером n есть x=na. Тогда
уравнение движения этого атома будет иметь вид
(8)
Решение
этого уравнения ищем в виде
(9)
Подстановка
этого выражения в уравнение приводит к необходимости соотношения
(10)
которое
является дисперсионным соотношением. Рис.2 представляет это соотношение.w, как видно, является периодической функцией k и период её
есть вектор обратной решетки 2p/a,
поэтому можно считать, что все состояния принадлежат первой зоне Бриллюэна
(выделена желтым на рисунке)
Рисунок
2 - Соотношение дисперсий.
При
k
0
,
т.е.
дисперсионное соотношение асимптотически переходит в таковое для сплошного
стержня.
Граничные
условия в форме Кармана - Борнаun+L=un, дают
eikLa =1, т.е. 
.
Соседние
значения k отличаются друг от друга на
.
Здесь
La - длина кристалла. В трехмерном случае это переходит
в
. Если через все значения ki провести плоскости, параллельные узловым плоскостям
обратной решетки, то первая зона Бриллюэна разобьется на N
фазовых ячеек с объёмом D3k, каждой из них принадлежит одно состояние. Т.о. в
одноатомном трехмерном кристалле осуществляется N продольных
колебательных состояний, кроме того, в таком кристалле осуществляется еще и 2N
поперечных колебаний
Рисунок 3 - Продольные и поперечные акустические колебания. Обозначения LA - продольная акустическая мода,TA - поперечная акустическая мода
1.4 Линейный кристалл с двумя атомами в
ячейке
Наличие двух атомов в ячейке в одномерном кристалле дает две ветви
колебаний. Пусть линейный кристалл с двумя атомами M1 и M2 в
элементарной ячейке подвержен возмущению упругими силами. Полагая силы
пропорциональными изменению расстояния между соседними атомами, а отклонения от
положений равновесия, обозначая u2n, u2n+1, u2n+2,...., запишем уравнения движения атомов с номерами 2n и 2n+1, имеющими массы M1иM2 ,соответственно:
(11)
Решение
будем искать в виде:
(12)
Подстановка
этих решений в уравнения приводит к системе алгебраических линейных уравнений
(13)
Эта
система имеет нетривиальное решение, когда определитель, составленный из
коэффициентов при неизвестных u1 иu2, равен нулю:
Откуда
получаем:
(14)
Если
перед радикалом взять “-“, то получаем т.н. акустическую ветвь колебаний (LA).
В этом случае при k=0 получим- w-=0 и u1=u2, а
при k= p/2a имеем w-=Ö2b/M2. Если взять “+”, то при k=0 w+=
, а -M1u1=M2u2, а
при k= p/2aw+=Ö2b/M1. Эта ветвь называется продольной оптической (LO),
так как в случае ионных кристаллов из-за противофазного движения атомов с
разной массой возникает переменный дипольный момент у элементарных ячеек, что
приводит к взаимодействию кристалла с электромагнитной волной.
Рисунок
4 - Оптическая и акустическая мода колебаний.
Обозначения
LO - продольная оптическая мода, LA -
продольная акустическая мода.
Налагая граничные условия Кармана-Борна, получаем значения
k=(2p/2aL)m,
где m=±1,±2,±3,.....,±L/2, где L - число элементарных ячеек в кристалле. Т.о. в данном
кристалле имеется 2Lколебательных
состояний. Если учесть и поперечные колебания, то добавляются две акустически и
две оптические ветви, точнее в нашем случае одна дважды вырожденная
акустическая ветвь и одна дважды вырожденная оптическая ветвь. В этом случае
общее число состояний 3*2*L.
Обобщая на трехмерный кристалл, получаем общее число ветвей 3*n, из которых три акустические, а
остальные - оптические (n-
число атомов в элементарной ячейке). Отсюда следует, что общее число
колебательных состояний в трехмерном кристалле равно числу степеней свободы, то
есть 3nL
Рисунок 5 - Структура колебательных мод сложного кристалла.
Колебания решетки проявляются в рэлеевском и комбинационном рассеянии
света, в поглощении и испускании ИК-излучения (остаточные лучи), их энергия
является основной частью внутренней энергии кристалла. Спектры колебаний
получают, исследуя форму дифракционных спектров рентгеновских лучей,
нейтронограмм и др.
Колебания решетки рассматривается как система осцилляторов.В том
приближении, которое мы рассмотрели, колебания решетки не зависят друг от друга
- это нормальные колебания, полностью определяемые величинами u0, k,w. Квантово механически это означает, что каждое колебание представляет
собой одномерный гармонический осциллятор с частотой w, возбужденный на энергетический
уровень v, т.е. имеющий энергию
(15)
причем w=w(l,k), l- номер ветви, а v=0,1,2,3,4,..... Т.о. энергия всех
колебаний есть
(16)
.5 Обзор теорий колебаний кристаллической решётки с дефектами
В настоящее время существует определенная
непоследовательность в подходах к описанию динамики решетки идеального ионного
кристалла и кристалла с дефектом. При анализе колебаний идеального кристалла
общепризнанным является учет поляризуемости ионов, что в рамках оболочечной
модели обеспечивается рассмотрением движения остовов и оболочек ионов. При
анализе локализованных колебаний кристалла, однако, рассматривается только
движение остовов ионов. Эта непоследовательность оправдана в тех случаях, когда
необходимо получить, лишь качественное представление о колебаниях решетки
вблизи дефекта и когда поляризуемость ионов невелика. Такой подход позволяет
сильно упростить задачу. Однако при количественном анализе экспериментальных
данных этот подход может оказаться неприемлемым, поскольку при изучении многих
физических процессов, таких как ИК-поглощение света, комбинационное рассеяние
света, электронно-колебательные переходы и т.д. интерес представляет поведение
внешних электронов, составляющих подвижные оболочки ионов.
Рассмотрим построение матрицы возмущения
для динамических уравнений оболочечной модели и вычислим функции Грина для
колебаний остовов и оболочек ионов примесного кристалла.
. Функции Грина
Для получения всех необходимых характеристик дефектных колебаний
достаточно вычислить соответствующие фононные функции Грина. Матрица функций
Грина дефектного кристалла 
вычисляется через матрицу аналогичных функций идеального кристалла
методом, развитым Лифшицем
(17)
Здесь 
-матрица возмущения динамической матрицы Lидеального кристалла при внесении дефекта
(18)
E-единичная матрица.
Уравнения динамики оболочечной модели в
стандартном виде выглядят следующим образом
(19)
где U - вектор смещений остовов, a W-вектор
смещений оболочек относительно остовов ионов; М - диагональная матрица масс
ионов. Матрицы силовых постоянных 
- выражаются через матрицы близкодействия Rи кулоновского взаимодействия С
(20)
где Z и Y-диагональные матрицы
зарядов остовов и оболочек ионов соответственно, k-диагональная матрица
констант связи остов-оболочка. Уравнение для смещений остовов ионов U.
(25)
где динамическая матрица 
имеет достаточно сложный вид
(22)
Для изовалентной примеси замещения происходит локальное возмущение
близкодействия в ее окрестности
(23)
т.е. матрица 
имеет конечный размер. Однако вследствие
сложной зависимости (6) динамической матрицы 
от матрицы близкодействия Rматрица 
становится бесконечной.
Чтобы обойти выше указанную трудность, оболочки ионов
рассматриваются как независимые составляющие кристаллической решетки с
маленькой массой т. Динамическая матрица такой задачи имеет вид
(24)
Если возмущение затрагивает п типов колебаний дефектной области
кристалла, то с учетом движения оболочек ионов количество колебаний
удваивается, а количество функций Грина, которое необходимо вычислить, равно 
.
Вторая трудность связана с решением задачи на собственные значения
с динамической матрицей (28), но высокочастотные колебания оболочек ионов
создают эффективную диэлектрическую проницаемость среды для остовов ионов. С
математической точки зрения это сводится к модификации только вещественных
частей функций Грина колебаний оболочек ионов, в которых появляются постоянные
добавки в фононной области частот неисчезающие при 
. Для их вычисления достаточно выполнить
очень грубое решение задачи с матрицей (28). В остальном вид функций Грина в
фононной области частот остается без изменений, поэтому подробные вычисления
частот и векторов поляризации нужно проводить лишь для фононной области
спектра. При этом для вычисления колебаний остовов ионов достаточно использовать
уравнение (25), а смещения оболочек находить из второго из уравнений (23).
Вычисления
функций Грина дефектных колебаний, как известно значительно упрощаются, если
перейти от удельных смещений к их линейным (симметризованным) комбинациям Q(Г),
преобразующимся по неприводимым представлениям Г точечной группы симметрии
примесного центра
(25)
где Т(Г)-прямоугольная матрица преобразования.
Полносимметричное возмущение решетки может смешивать только
колебания одинаковой симметрии, поэтому выражение (1) распадается на
независимые выражения для каждого представления Г
(26)
где 
-матрица функций Грина
смещений симметрии Г. а дефектная матрица 
вычисляется по формуле
(27)
Размерность матрицы 
определяется количеством смещений одинаковой симметрии Г в
дефектной области кристалла и при высокой симметрии позиции примеси значительно
ниже размерности матрицы 
. Из общего матричного выражения (30)
легко получить пригодное дляпрактических расчетов выражение для функций Грина 
смещений i и j
(28)
Если же одно из смещений (например, i) или оба принадлежат к
дефектной области кристалла, то выражение упрощается
(29)
Методы расчета функций Грина смещений идеального кристалла,
входящих в выражения (32), (33), хорошо известны.
Вещественные части функций Грина
вычисляются из дисперсионных соотношений.
Для высокочастотных колебаний, необходимо
вычислить лишь вещественную часть функций Грина смещений оболочек при нулевом
значении частоты
(30)
Суммирование в (34) проводится только по высокочастотным модам
колебаний.
Для функций Грина смещений остовов ионов существуют простые
условия нормировки
(31)
Эти условия выполняются и для функций
Грина дефектных колебаний, что дает хороший способ контроля правильности
вычислений.
Для смещений оболочек ионов таких условий
нормировки нет, поскольку основная интенсивность их функций Грина сосредоточена
в высокочастотной области спектра.
Выражения (19)-(35) достаточны для решения поставленной задачи,
если определена модель дефекта 
.
2. Модель дефекта
Изовалентное замещение основного иона
решетки примесным ионом может сопровождаться следующими статическими и
динамическими искажениями решетки.
а)Полносимметричной деформацией решетки,
затрагивающей одну или несколько координатных сфер примеси.
б)Статической поляризацией ионов
ближайшего окружения (лигандов) примеси.
в)Изменением короткодействующих сил
взаимодействия между примесным ионом и лигандами, а также между основными ионами
решетки в окрестности примеси.
г)Изменением динамических свойств решетки,
связанным с отличием поляризуемости примесного иона от поляризуемости
замещенного иона решетки.
д)Изменением динамических свойств решетки,
связанным с отличием массы примесного иона от массы замещенного иона решетки.
Здесь мы рассмотрим влияние всех перечисленных факторов на
динамику решетки, за исключением статической поляризации лигандов. Ее учет не
представляет принципиальных затруднений, однако требует перехода к другим динамическим
уравнениям с абсолютными смещениями остовов 
и оболочек 
ионов. Полносимметричная деформация
решетки обусловливает изменение кулоновского взаимодействия 
, а также часть изменения близкодействия 
за счет изменения расстояний между
ионами. Другая часть 
обусловлена отличием электронной оболочки
примесного иона от оболочки замещаемого (основного) иона. Различие
поляризуемостей примесного и основного ионов описывается изменением заряда
оболочки 
и константы связи остов-оболочка 
. Однако, поскольку для оболочечной модели
характерна тесная корреляция между Yи k, можно рассматривать только изменение . Наконец, различие масс примесного и
основного ионов дает дефект массы 
. В рамках рассмотренной модели легко, исходя из общего вида
динамических уравнений (3), написать выражение для дефектной матрицы
(32)
Дальнейшая детализация выражения для возможна для конкретного типа кристалла.
1.6 Теория зон Бриллюэна
Бриллюэна
зоны, области значений волнового вектора
<#"705714.files/image001.gif"> определяется формулой
Где
h-элемент
симметрии фактор-группы кристалла;
- число атомов, которые с точностью до вектора примитивной
трансляции не изменяют своего положения в результате действия операции
симметрии h;
- угол поворота при действии операции симметрии h.
При этом знак плюс соответствует собственным вращениям, знак минус
- несобственным вращениям.
Характер полного колебательного представления 
содержит характеры акустических 
и оптических 
колебаний:
,
(34)
С помощью формулы приведения характер представления оптических
колебаний можно разложить по неприводимым представлениям точечной группы 
кристалла:
, 
,(35)
где mi- число, которое показывает, сколько раз i-е неприводимое представление группы кристалла содержится в колебательном представлении число элементов симметрии (порядок) в группе ; 
- характер - j-го неприводимого представителя группы для операции симметрии hj; Nj-
число операции j-го класса симметрии.
Различные
примеры классификации оптических колебаний в кристаллах с использованием общего
метода в работах. Наиболее трудным моментом в этом методе является определение
числа атомов , не изменяющих своего положения при действии операции
симметрии hфактор-группы кристалла. Метод допускает проведение
классификации колебаний по типам симметрии при неизвестной позиционной
симметрии атомов, однако если позиционная симметрия известна, то решение этой
задачи существенно упрощаются.
Для
теоретико-групповой классификации колебаний решетки кристаллов, содержащих в
структуре многоатомные изолированные группы (молекулы), более удобен метод
позиционной симметрии. В этом методе анализируется влияние кристаллического
поля конкретной симметрии на нормальные колебания молекулы (комплексного иона),
находящийся в свободном состоянии, путем корреляции колебательных уровней
кристалла с колебательными уровнями свободной молекулы. Схема анализа сводится
к следующему. Симметрия молекулы в свободном состоянии характеризуется
некоторой точечной группой симметрии. В кристаллической структуре молекула, как
правило, занимает позицию с более низкой симметрией. Эта особенность может
привести к тому, что в кристалле под действием кристаллического поля может быть
снято вырождение внутренних колебаний, имеющее место для молекулы в свободном
состоянии. При этом может возникнуть расщепление некоторых колебаний. Такое
расщепление колебаний получило название Бете-расщепления. Более сильное
воздействие на молекулы, но и сдвиг частот, а также изменение правил отбора. В
результате в спектре кристалла могут проявиться некоторые типы внутренних
колебаний, запрещенные для свободной молекулы.
Для
учета действия кристаллического поля на колебательные уровни внутренних
колебаний следует колебания кристалла классифицировать по неприводимым
представлениям группы позиционной симметрии. На практике это достигается с
помощью корреляции между неприводимыми представлениями точечной группы
симметрии свободной молекулы и группы ее позиционной симметрии в элементарной ячейке
кристалла. Если элементарная ячейка кристалла содержит не одну молекулу, а
несколько, то кроме действия статического кристаллического поля необходимо
учитывать динамическое взаимодействие одинаковых молекул в ячейке, которое
приводит к давыдовскому расщеплению внутренних колебаний молекулы в кристалле.
Число компонент давыдовского расщепления равно числу молекул в примитивной
ячейке. Величина давыдовского расщепления является мерой межмолекулярного
взаимодействия. Это расщепление наиболее значительно в молекулярных и
молекулярно-ионных кристаллах, где оно составляет обычно 10-15см-1 .
Первоочередной
задачей при классификации колебаний в методе позиционной симметрии является
выяснение колебательного представления молекулы в свободном состоянии.
Колебательное движение изолированных молекул различной симметрии в свободном
состоянии в настоящее время в литературе подробно изучены. Молекула, состоящая
из z атомов, имеет 3z колебательных
степеней свободы. Из них три колебания соответствуют трансляционному движению
молекулы как целого, три колебания - вращениям молекулы как целого. При этом
для линейных молекул имеется только две вращательные степени свободы. Остальные
3z-6 колебаний являются внутренними колебаниями атомов в
молекуле.
Корреляции
между неприводимыми представлениями точечной группы симметрии молекулы в
свободном состоянии и группы ее позиционной симметрии в кристаллической
структуре осуществляются с помощью таблиц корреляций типов симметрии групп и
подгрупп. Подобные таблицы и схемы проведения анализа для конкретных кристаллов
приведены во многих монографиях и обзорных статья. Суть проведения корреляций
сводится к следующему. Если точечная группа симметрии Q является
подгруппой другой точечной группы P, может быть разложено на
неприводимые представления группы Q по формуле, где f -
число элементов (порядок) группы Q. Формула устанавливает корреляцию i-го
неприводимого представления группы P с j-ми
неприводимыми представлениями подгруппы Q. Эта
корреляция устанавливается с помощью таблицы характеров неприводимых
представлений точечных групп симметрии. Некоторые точечные группы симметрии
кристалла имеют несколько подгрупп с одинаковым набором элементов симметрии, но
различающиеся ориентацией элементов. В таких ситуациях возникает
неоднозначность корреляции неприводимых представлений групп и подгрупп. Эта
неоднозначность может быть устранена, если установлена (например,
дифракционными методами анализа структуры) точная ориентация молекулярной
группы в элементарной ячейке кристалла.
2.
МЕТОДЫ РАСЧЕТА ЭНЕРГИИ ВЗАИМОДЕЙСТВИЯ АТОМОВ В КОВАЛЕНТНЫХ КРИСТАЛЛАХ
2.1 Расчет энергии взаимодействия из геометрических соображений
Наиболее простой способ расчета энергии межатомного взаимодействия
заключается в использовании следующей модели. Предполагается, что атомы,
находящиеся в узлах кристаллической решетки представляют собой положительно
заряженные точечные ядра, окруженные отрицательно заряженным облаком
электронов, плотность которых считается постоянной. Радиус электронного облака
полагаем равным радиусу иона. Такая модель представляет собой довольно грубое
упрощение модели Томаса-Ферми. Тем не менее, такая модель представляет
определенный интерес.
Найдем
вначале размеры области перекрытия двух однородно заряженных шаров, имеющих
радиусы 
и 
.
Расстояние между центрами шаров обозначим через 
. Тогда
(см. рис.6):
Рисунок
6
Условием
перекрытия шаров является неравенство 
. Из
чертежа находим: 
, 
, откуда 
. С другой стороны, 
.
Следовательно, 
, тогда расстояния 
, 
можно найти из системы уравнений:
(36),
(37)
Отсюда
, 
.
Выражение для 
можно записать в виде:
(38)
Последнее
выражение также можно записать в виде:
или
(39)
Площадь
перекрытия шаров в сечении r равна:
(40)
В
работах принято, что энергия взаимодействия пропорциональна площади. Оценим эти
величину площади S для атомов кристаллической решетки кристалла
вольфрамата бария BaWO4.
На
наш взгляд, более последовательный подход к вычислению энергии взаимодействия
состоит в том, что она полагается равной суммарной энергии кулоновского
взаимодействия атомов:
(41)
где
первое слагаемое отвечает энергии взаимодействия двух точечных ядер с зарядами
и 
, второе
слагаемое отвечает взаимодействию электронной оболочки, несущей суммарный заряд
с ядром, заряд которого равен (в силу сферической симметрии электронная оболочка
первого атома создает в месте расположения ядра второго атома потенциал 
). Третье слагаемое отвечает взаимодействию
электронной оболочки второго атома с ядром первого атома. Энергия 
, связанная с перекрытием электронных оболочек, равна:
( 42)
где
, 
, 
, 
, модуль
радиус-вектора 
изменяется в пределах 
, угол 
изменяется в пределах 
,
величины 
, 
изменяются
в аналогичных пределах. Расстояние 
равно
(43)
2.2 Расчет энергии взаимодействия
методом Томаса-Ферми-Дирака-Гамбоша
Наиболее точными количественными расчетами
энергии взаимодействия атомов гетеродесмических соединений в кристаллической
решетке являются метод Хартри-Фока и метод функционала плотности, базирующийся
на уравнениях Кона-Шема. Эти методы освещены в обзорах Н. Марча, В. Кона, П.
Вашишта, С. Лундквиста, А. Уильямса, У. Барта, Н. Ланга. К недостаткам этих
методов следует отнести значительную трудоемкость расчетов, которые, кроме
того, являются «непрозрачными»: правильность расчетов, как правило, можно
установить только после того, как расчеты выполнены. Действительно, метод
Хартри-Фока базируется на численном решении нелинейных уравнений, и результаты
промежуточных вычислений дают мало информации об окончательном результате.
Наряду с перечисленными методами существует значительно более
простой метод, основанный на статистической теории атома
Томаса-Ферми-Дирака-Гамбоша. Расчеты в рамках этого метода базируются на
сравнительно простых уравнениях, допускающих ясную физическую интерпретацию на
всех этапах вычислений. В рамках этого метода энергия связи между двумя атомами
складывается из кулоновского взаимодействия остовов атома 
, энергии электростатического взаимодействия
электронов ие, корреляционной энергии вследствие взаимодействия электронов с
антипараллельными спинами 
, энергии экранирования ионов 
, кинетической uk и обменной ua энергий электронов
в объеме межатомной связи:
(44)
Энергия кулоновского взаимодействия обращается в нуль при
взаимодействии нейтральных атомов одного сорта, для атомов разного сорта эта
энергия не равна нулю. Сумма 
по модулю практически равна половине кинетической энергии
электронов, образующих связь 
. В отличие от положительной кинетической энергии обменная энергия
ua и энергия us имеют отрицательный знак.
Обменная энергия ua и
кинетическая энергия uk электронного
газа выражаются через электронную плотность 
в объеме связи. Как было показано в, градиентные поправки
неоднородного электронного облака связи малы, и ими можно пренебречь.
В приближении Томаса-Ферми-Дирака-Гамбоша парная связь атомов A-B,
обусловленная кинетической uk, обменной ua и корреляционной 
энергиями электронного газа в объеме
связи, равна
(45)
где кк = 2.871 а т.ед., ка = 0.7381 ат.ед. (без учета 
), и ка = 0.8349 ат.ед. (с учетом ).
На половине длины парной связи 
атомов A, B в точках, равноудаленных от центров этих атомов, 
.
Объемная плотность суммарной энергии в данной точке на половине
длины связи:
(46)
Результирующие электронные плотности орбиталей, образующих связь,
описываются либо в форме выражений Слейтера, либо в форме функции Гаусса,
пользующейся широким распространением. В этом случае электронная плотность в
зависимости от расстояния r от центра атома
может быть представлена выражением 
(коэффициент 
пропорционален обратной величине квадрата радиуса связывающей
орбитали). При N - числе электронов атома, участвующих в
связи, предэкспоненциальный множитель А определяется из условия нормировки
(47)
На половине длины связи расстояние r до любой точки электронного облака связи, расположенной на
нормали к связи на расстоянии у от нее, определяется из условия 
.
Тогда линейная плотность энергии связи на половине ее длины 
равна:
(48)
Распределение линейной плотности энергии связи вдоль линии связи l аппроксимируется функцией вида 
, в которой 
. При этом энергия парной связи 
определяется интегрированием линейной плотности вдоль линии связи
(49)
после введения коэффициента 
часто учитывающего роль 
полная энергия парной связи в данной координатной сфере i кристалла элемента j будет
(50)
и при полном учете корреляционной энергии 
(51)
Здесь 
- длина связи между данным атомом с
атомами первой координационной сферы, Ki -
отношение радиуса iкоординационной сферы к радиусу первой
сферы r1, 
элемента j. По
существу, численные коэффициенты в выражении 
являются вариационными параметрами 
, 
(52)
Каждый атом элемента jвзаимодействует
с 
атомами i координационной сферы. При этом
меньше 
соответственно коэффициенту экранирования
(53)
- относительное эффективное сечение атома (
).
Полная энергия связи кристалла, приходящаяся на один атом элемента
j, равна
(54)
3.
РАСЧЁТ СПЕКТРА КОЛЕБАНИЙ КРИСТАЛЛИЧЕСКОЙ РЕШЁТКИ ВОЛЬФРАМАТА БАРИЯ
.1 Кристаллическая структура вольфрамата бария
Вольфраматы и молибдаты щелочноземельных металлов являются перспективными
материалами квантовой электроники используемыми в качестве преобразователей
частоты лазерного излучения на основе эффекта вынужденного комбинационного
рассеяния (ВКР), т.е. соединения, имеющие формулу BaWO4, Эти кристаллы относятся к классу соединений с
гетеродесмической (смешанной ионно-ковалентной) структурой, что определяет
особенности их фононного спектра и, как следствие, использование их в качестве
материалов для ВКР- преобразователей.
Монокристаллы вольфрамата бария относятся к тетрагональной сингонии. Их
симметрия описывается пространственной группой I 41/a. Параметры элементарной
ячейки BaWO4: а=5.614 Å, b=5.614 Å, с=12.719 Å, б=в=г=90 град, объем элементарной
ячейки 400.86 Å3 [Быличкина Т.И., Солева Л.И., Победимская Е.А. и др.
Кристаллическая структура Ва молибдата и Ва - вольфрамата // Кристаллография.
1970. Т. 15. № 1.С. 165-167.].
Координаты
базисных атомов элементарной ячейки представлены в таблице 3. Зона Бриллюэна
для кристалла BaWO4 имеет симметрию той же группы что и сам кристалл, и
характеризуется следующими базисными векторами: 
,
, 
.
Для
расчета фононных спектров вольфрамата бария использовался метод, описанный в
работе [4]. Потенциал межатомного взаимодействия 
выбирался
в виде:
(55)
суммирование распространяется на все пары ионов, принадлежащих
элементарной ячейке. Потенциальная энергия 
взаимодействия ионов б и в,
находящихся на расстоянии
друг от друга, выбиралась в виде:
(56)
где
,
- векторы смещений атомов б и в из положений
равновесия, определяемых векторами 
, 
; 
; 
- равновесное значение расстояния между атомами б и
в, 
- обрезающий множитель, учитывающий короткодействие
потенциала. При расчете фононного спектра обрезание потенциала удобнее
проводить в k-пространстве волновых векторов.
Поскольку кристаллы вольфраматов имеют смешанный тип связи, то потенциал
взаимодействия 
выбирался в виде потенциала Морзе
(57)
Потенциал Морзе и все его производные регулярны при любых расстояниях
между ионами, что выгодно отличает его от сингулярного дальнодействующего
кулоновского потенциала и его модификаций, учитывающих экранирование.
Характерное для потенциала Морзе быстрое убывание энергии взаимодействия с
возрастанием расстояния между ионами позволяет достаточно адекватно описывать
гетеродесмические кристаллы, у которых имеет место чередование сильных ионных и
слабых ковалентных связей. Обоснование потенциала Морзе можно получить в рамках
теории Томаса-Ферми-Гамбоша. Параметр Ω, входящий в (57), является
подгоночным параметром, и его величина определяется по краю инфракрасной полосы
поглощения.
Таблица №3. Координаты базисных атомов элементарной ячейки BaWO4
|
Атом
|
Координаты базисных атомов
|
|
x/a
|
y/b
|
z/c
|
|
Ba
|
0.5
|
0.75
|
0.125
|
|
W
|
0.5
|
0.25
|
0.375
|
|
O(1)
|
0.646
|
0.146
|
0.291
|
|
O(2)
|
0.397
|
0.104
|
0.459
|
|
O(3)
|
0.604
|
0.396
|
0.459
|
|
O(4)
|
0.354
|
0.354
|
0.291
|
Матрица
силовых постоянных 
, вычисленная в декартовых координатах x, y, z
вектора смещения равна:
. (58)
Фононный
спектр определялся путем решения секулярного уравнения для Фурье-компонента 
вектора смещения:
(59)
где
динамическая матрица D связана с матрицей силовых постоянных 
равенством
.
При
расчете фононного спектра BaWO4 величина параметра Ω, определенная по
ИК-границе полосы пропускания, взята равной 21.6. Параметр для каждой пары ионов выбран как равновесное
расстояние между ними. Такой выбор соответствует минимуму потенциальной энергии
(
) и
приводит к значительному упрощению расчетов. Параметры межатомного
взаимодействия 
выбраны пропорциональными суммарным значениям соответствующих
потенциалов ионизации для соответствующих пар ионов, их численные значения
приведены в таблице 3.
Таблица
4. Параметры межатомного взаимодействия 
для пар
ионов BaWO4.
|
Пара ионов
|
Параметрмежатомн. взаимодействия (эВ)Пара ионовПараметрмежатомн. взаимодействия (эВ)
|
|
|
|
Ba-W
|
3.5
|
W-O
|
5.1
|
|
Ba-O
|
4.7
|
O-O
|
4.7
|
Для
определения законов дисперсии 
секулярное
уравнение (59) записывалось в k-пространстве при помощи косинус - преобразования
Фурье. Короткодействующий характер взаимодействия учитывался следующим образом:
при нахождении Фурье-образа, матрицы силовых постоянных, для которых 
, полагались равными нулю. Тем самым предполагалось,
что потенциал взаимодействия простирается на расстояние меньшее половины
размера элементарной ячейки. По этой же причине, при вычислении Фурье-образа
матрицы силовых постоянных, Фурье-преобразование проводилось только в пределах
элементарной ячейки.
Диагонализация
силовой матрицы в k-пространстве производилась численно с помощью
программы Maple 4.0 для сечения 
, в
котором было выбрано 6 точек с шагом 
.
Сходимость численного расчета обеспечивалась путем добавления атомов в
элементарную ячейку, (последовательным увеличением размера силовой матрицы), и
сравнением влияния добавления на вычисленные ранее собственные значения. Расчет
считался оконченным, когда относительное изменение собственных значений
составляло величину порядка 10-6. Реальная точность расчета спектра колебаний
определяется, во-первых, аппроксимацией потенциала межатомного взаимодействия,
и, во-вторых, обрезанием матричных элементов в k-пространстве.
Влияние описанных факторов снижения точности наиболее сильно в низкочастотной
области спектра, где погрешность в энергии фонона составляет величину порядка
самой энергии.
Диэлектрические
кристаллы, активированные ионами лантаноидов и переходных металлов представляют
большой интерес в квантовой электронике. К таким материалам относятся, среди
прочих, кристаллы вольфраматов, активированные ионами Cr3+.
Представляет интерес расчет фононного спектра кристаллов BaWO4:Cr3+,
который предполагается выполнить в дальнейшем.
.2 Описание программы для расчёта частоты колебаний кристаллической
решётки
Для расчёта частоты колебаний кристаллической решётки взяли за центр атом
скандия и с помощью программы Topos,
определили ближайшее окружение кислорода и вольфрама, в котором выделили
эквивалентные позиции, отличающиеся на трансляционный период и из остальных
позиций взяли количество соответствующих химической формуле 1-Ba,2-W, 3,4,5,6-O, и для этих атомов определили их
координаты в ячейке. Взяли потенциал Морзе, который предлагается в статье
Дубровского и на его основе рассчитали силовые константы, составили из них
силовую матрицу, для которой составили секулярное уравнение для определения
частоты колебаний от центра зоны Бриллюэна от 0 через 0,1 до 0,5.
Рисунок 7- Элементарная ячейка кристалла вольфрамата бария.
Рисунок 8 - Фрагмент кристаллической структуры BaWO4, использованный для расчета фононного спектра
(62)
где
суммирование распространяется на пары атомов, занимающих неэквивалентные
позиции, потенциальная энергия 
равна
(63)
-
расстояние между атомами 
и 
, которое
вычисляется по данным таблицы 1,
-
параметр межатомного взаимодействия между атомами и ,
пропорциональный энергии взаимодействия атомов, приведенный в таблице 2(общий
множитель пропорциональности включается в определение 
),
, 
- равновесное значение расстояния,
-
подгоночный параметр, величина которого определяется по длинноволновому краю
полосы прозрачности кристалла BaWO4.
Вариация
потенциала межатомного взаимодействия 
позволила
получить матричные элементы динамической матрицы
(здесь
- массы атомов , , 
-
декартовы компоненты вектора смещения -го атома
из положения равновесия 
), являющиеся основой для составления секулярного
уравнения для Фурье-компонентов вектора
(64)
для
численного решения которого была составлена программа на языке программирования
Марle.
Рисунок
9 - Кривые дисперсии, содержащие одну троекратно вырожденную акустическую ветвь
и пять троекратно вырожденных оптических ветвей. Энергия колебаний
(вертикальная ось) рассчитана в ридбергах Ry: 1 Ry=13.6
эВ. По горизонтальной оси отложена величина волнового числа вдоль направления 
в долях 
.3 Расчет плотности фононных состояний для вольфрамата бария
Приведенный на рис.9 фононный спектр BaWO4 использовался для расчета плотности фононных
состояний численным методом по формуле
,
где
- объем элементарной ячейки. Диапазон энергий
0...0,45 эВ был разбит на 550 интервалов, и для каждого интервала была найдена
длина соответствующего интервала на отрезке 
. Эта
длина пропорциональна искомой плотности состояний g(e). Результаты расчета плотности состояний представлены на рисунке 10.
Рисунок
10 - Гистограмма, результаты расчета плотности фононных состояний.
В представленных на рис.10 результатах расчета плотности фононных
состояний обращает на себя внимание, что плотность числа колебательных
состояний значительно возрастает для первой(акустической) и пятой(оптической)
ветвей. Нумерация ветвей производится снизу вверх.
Заключение
Расчеты в рамках предложенной модели показывают что спектр
делокализованных колебаний кристаллической решетки вольфрамата бария состоит из
отдельных изолированных пиков, т.е. ведет себя подобно спектру локализованных
колебаний изолированной молекулы. Такая локализация колебаний вполне объясняет
эффективность вольфрамата бария как преобразователя частоты лазерного излучения
за счет эффекта комбинационного рассеяния, поскольку известен факт высокой
эффективности молекулярных сред как ВКР-преобразователей.
Список использованных источников
1. Lebedev V.A., Pisarenko V.F., Selina N.V.,
Perfilin A.A., Brik M.G. Spectroscopic and luminescent properties of Yb, Er:
LaSc3(BO3)4 crystals. //Optical Materials. 2000 14 C.121-126.
. Selina N.V., Avanesov A.G., Lebedev V.A.,
Stroganova E.V., Tumayev E.N., Brik M.G., Special features of the phonon
spectrum and non-radiative transitions in the Cr3+-doped ionic-covalent
crystals. // Solid State Communications 2008, 146, 298.
3. ГорюновА.В., Кузьмичева Г.М.Рентгенографическое
исследование кристаллов LaSc3(BО3)4, активированных ионами хрома и неодима //
Журнал неорганической химии. 1996. Том 41, №10. С.1605-1611.
. Дубовский О.А., Орлов А.В., Семенов В.А.
Колебательные спектры металл-металлоидных кристаллитов Fe3P, Fe2P: фононные и
бризерные возбуждения. // ФТТ. 2003. Том 45, Вып.2. С.309-316.
. Ландау Л.Д., Лифшиц Е.М. Квантовая
механика:нерелятивистская теория. М.:Наука, 1974,752 с.
. Гамбош П.Статистическая теория атома и ее
применения. М:, Иностранная литература, 1951. 398 с.
. Рыбаков В.Б., Кузьмичева Г.М., Мухин Б.В., Жариков
Е.В., Агеев А.Ю., Кутовой С.А., Кузьмин О.В. Рентгеноструктурное исследование
соединений (Ce,Gd)Sc3(BO3)4 семейства хантита // Журнал неорганической химии.
1997. Том 42, №1. С.9-16.
. Грибов.Л.А., Баранов В.И., Зеленцов Д.Ю.Электронно-колебательные
спектры многоатомных молекул. М.: Наука, 1997. С.197-201.
. Акопян И.Г., Писаренко В.Ф., Тумаев Е.Н. Фононный
спектр кристаллов скандобората лантана (LSB) и кристаллов LSB, активированных
ионами Cr3+ // Экологический вестник научных центров Черноморского
экономического сотрудничества, 2005, №3, С.60-65.
. Строганова Е.В. Спектрально-люминесцентные и
кинетические исследования кристаллов двойных церий-скандиевых боратов
(Ce,Gd)Sc3(BO3)4, активированных ионами Cr3+: Дис. канд. физ.-мат. наук, Краснодар,
2003. 102 с
. Пуле А., Матье Ж.-П. Колебательные спектры и
симметрия кристаллов.-М.: Мир, 1973.-437 с.
. Рассеяние света в твердых телах / Под ред.
М.Кордоны.-М.: Мир, 1979.-379 с.
. Колебательная спектроскопия. Современные воззрения.
Тенденции развития / Под ред. А. Барнса и У. Орвилл-Томаса. - М.: Мир, 1981. -
480 с.
. Сущинский М.М. Спектры комбанационного рассеяния
света молекул и кристаллов. - М.: Наука, 1969.-576 с.
. Жижин Г.Н., Маврин Б.Н., Шабанов В.Ф. Оптические
колебательные спектры кристаллов. - М.: Наука, 1984.- 232 с.
. Боборов А.В., Мулдахметов З.М. Спектроскопия
комбинационного рассеяния света. - Алма-Ата: Наука, 1981. - 152 с.
. Горелик В.С., Умаров Б.С. Введение в спектроскопию
комбинационного рассеяния света. - Душамбе: Дониш, 1982. - 287 с.
. Китайгородский А.И. Молекулярные кристаллы. - М.:
Наука, 1971.- 424 с.
. Применение спектров комбинационного рассеяния / Под
ред. А.Андерсона и К.И. Петрова. - М.: Мир, 1977. - 586 с.
. Борн М., ХуангКунь. Динамическая теория
кристаллических решеток.- М.: ИЛ, 1958. - 488 с.
21. Zhizhin G.N., Muktarov E.I. Optical spectra
and lattice dynamics of molecular crystals. Vibrational spectra and structure /
Eds. J.R. Durig. A series of advances, V. 21.- Amsterdam: ELSEVIER, 1995.- 490
p.
22. ЛисицаМ.П., ЯремкоА.М. РезонансФерми.- Киев: Науковадумка, 1984. - 334 с.
23. Шаскольская М.П. Кристаллография. - М.: Высшая школа.
1978.- 391с.
24. Horning D.F. // J.Chem.Phys. - 1948. - V.16,
№11. - P. 1063-2076.
25. Теория
неоднородного электронного газа / Под ред. С. Лундквиста и Н. Марча. Мир, М.
(1987). 400 с.
26. П. Гамбош.
Статистическая теория атома и ее применение / Под ред. К.П. Гурова. Изд-во
иностр. лит., М. (1951). 399 с.
. И.К. Дмитриева,
Г.И. Плиндов. Свойства атомов и ионов в свете статистической теории. НиТ, Минск
(1991). 189 с.
28. N.N. Sirota, I.M.
Sirota. Cryst. Res. Techn. 32, 1, 143 (1997).
29. В. Кон, П.
Вашишта. В сб.: Теория неоднородного электронного газа. Мир, М. (1987). С. 86.
. Дж. Эмсли.
Элементы. Мир, М. (1993). С. 256.
. М. Борн, Х.
Кунь. Динамическая теория кристаллических решеток. ИЛ. М. (1958). 488 с.
32. А.Марадудин,.Э.Монтролл,Д.Вейсс.
Динамическая теория кристаллических решеток в гармоническом приближении. Мир,
М. (1965). 384 с.
33. N.N. Sirota, T.D.
Sokolovcky. Chemical bond in Semiconductor and Solids. Consultants Bureau, N.Y.
(1967). P. 139.
34. S.C. Upadhyaya, J.C. Upadhyaya, R. Shyam. Phys. Rev.
B44, 1, 122 (1991).
35. Н.Н. Сирота, И.М.
Сирота, Т.М. Сошнина, Т.Д. Соколовский. ДАН 6 (2000).
36. Landolt-Bornsein.
Tabellen. Gruppe III B. Bd 13. Springer,N.Y. (1981).
37. H. Zeise.
Thermodynamik. Tabellen. S. Hirzel. Leipzig (1954).V. 3/1.
38. Физические
величины / Под ред. И.С. Григорьева, Е.З. Мейлихова. Энергоиздат, М. (1991).
Приложение
Программа для расчетов колебаний в кристаллической решётке вольфрамата
бария.
restart; a:=5.614; b:=5.614; c:=12.719;
> PP:=3.14159265358979; b1x:=2*PP/a; b1y:=0; b1z:=0;
b2x:=0; b2y:=2*PP/a; b2z:=0; b3x:=0; b3y:=0; b3z:=2*PP/a;
> Nfact:=0.0; OM:=21.6; y:=1.;
f1:=(2*OM+1)*exp(-2*OM*(y-1))-(OM+1)*exp(-OM*(y-1));
f2:=exp(-OM*(y-1))-exp(-2*OM*(y-1));
> qx:=Nfact*b1x; qy:=Nfact*b1y; qz:=Nfact*b1z;
> x1:=0.5*a; y1:=0.75*b; z1:=0.125*c; x2:=0.5*a;
y2:=0.25*b; z2:=0.375*c; x3:=0.646*a; y3:=0.146*b; z3:=0.291*c; x4:=0.397*a;
y4:=0.104*b; z4:=0.459*c; x5:=0.604*a; y5:=0.396*b; z5:=0.459*c; x6:=0.354*a;
y6:=0.354*b; z6:=0.291*c;
> r12:=sqrt((x2-x1)^2+(y2-y1)^2+(z2-z1)^2);
r13:=sqrt((x3-x1)^2+(y3-y1)^2+(z3-z1)^2);
r14:=sqrt((x4-x1)^2+(y4-y1)^2+(z4-z1)^2); r15:=sqrt((x5-x1)^2+(y5-y1)^2+(z5-z1)^2);
r16:=sqrt((x6-x1)^2+(y6-y1)^2+(z6-z1)^2);
> r34:=sqrt((x4-x3)^2+(y4-y3)^2+(z4-z3)^2);
r35:=sqrt((x5-x3)^2+(y5-y3)^2+(z5-z3)^2);
r36:=sqrt((x6-x3)^2+(y6-y3)^2+(z6-z3)^2);
> r45:=sqrt((x5-x4)^2+(y5-y4)^2+(z5-z4)^2);
r46:=sqrt((x6-x4)^2+(y6-y4)^2+(z6-z4)^2);
> r56:=sqrt((x6-x5)^2+(y6-y5)^2+(z6-z5)^2);
> gxx12:=(x1-x2)^2/r12^2; gyy12:=(y1-y2)^2/r12^2;
gzz12:=(z1-z2)^2/r12^2; gxx13:=(x1-x3)^2/r13^2; gyy13:=(y1-y3)^2/r13^2;
gzz13:=(z1-z3)^2/r13^2; gxx14:=(x1-x4)^2/r14^2; gyy14:=(y1-y4)^2/r14^2;
gzz14:=(z1-z4)^2/r14^2; gxx15:=(x1-x5)^2/r15^2; gyy15:=(y1-y5)^2/r15^2;
gzz15:=(z1-z5)^2/r15^2; gxx16:=(x1-x6)^2/r16^2; gyy16:=(y1-y6)^2/r16^2;
gzz16:=(z1-z6)^2/r16^2;
> gxx23:=(x2-x3)^2/r23^2; gyy23:=(y2-y3)^2/r23^2;
gzz23:=(z2-z3)^2/r23^2; gxx24:=(x2-x4)^2/r24^2; gyy24:=(y2-y4)^2/r24^2;
gzz24:=(z2-z4)^2/r24^2; gxx25:=(x2-x5)^2/r25^2; gyy25:=(y2-y5)^2/r25^2;
gzz25:=(z2-z5)^2/r25^2; gxx26:=(x2-x6)^2/r26^2; gyy26:=(y2-y6)^2/r26^2;
gzz26:=(z2-z6)^2/r26^2;
> gxx34:=(x3-x4)^2/r34^2; gyy34:=(y3-y4)^2/r34^2;
gzz34:=(z3-z4)^2/r34^2; gxx35:=(x3-x5)^2/r35^2; gyy35:=(y3-y5)^2/r35^2;
gzz35:=(z3-z5)^2/r35^2; gxx36:=(x3-x6)^2/r36^2; gyy36:=(y3-y6)^2/r36^2;
gzz36:=(z3-z6)^2/r36^2;
> gxx45:=(x4-x5)^2/r45^2; gyy45:=(y4-y5)^2/r45^2;
gzz45:=(z4-z5)^2/r45^2; gxx46:=(x4-x6)^2/r46^2; gyy46:=(y4-y6)^2/r46^2;
gzz46:=(z4-z6)^2/r46^2;
> e12:=3.5/13.6; e13:=e12; e14:=e12; e15:=3.5/13.6;
e16:=e15; e23:=5.1/13.6; e24:=5.1/13.6; e25:=4.7/13.6; e26:=4.7/13.6;
e34:=5.1/13.6; e35:=4.7/13.6; e36:=4.7/13.6; e45:=4.7/13.6; e46:=e45;
e56:=4.7/13.6;
> M1:=137.0; M2:=183.0; M3:=16.0; M4:=16.0; M5:=16.0;
M6:=16.0;> Dxx12:=2*OM*e12*(gxx12*f1+f2)/sqrt(M1*M2);
Dxx13:=2*OM*e13*(gxx13*f1+f2)/sqrt(M1*M3);
Dxx14:=2*OM*e14*(gxx14*f1+f2)/sqrt(M1*M4);
Dxx15:=2*OM*e16*(gxx12*f1+f2)/sqrt(M1*M5);
Dxx16:=2*OM*e16*(gxx12*f1+f2)/sqrt(M1*M6);
> Dxx23:=2*OM*e23*(gxx23*f1+f2)/sqrt(M2*M3);
Dxx24:=2*OM*e24*(gxx24*f1+f2)/sqrt(M2*M4);
Dxx25:=2*OM*e25*(gxx25*f1+f2)/sqrt(M2*M5);
Dxx26:=2*OM*e26*(gxx26*f1+f2)/sqrt(M2*M6);
> Dxx34:=2*OM*e34*(gxx34*f1+f2)/sqrt(M3*M4);
Dxx35:=2*OM*e35*(gxx35*f1+f2)/sqrt(M3*M5); Dxx36:=2*OM*e36*(gxx36*f1+f2)/sqrt(M3*M6);
> Dxx45:=2*OM*e45*(gxx45*f1+f2)/sqrt(M4*M5);
Dxx46:=2*OM*e46*(gxx46*f1+f2)/sqrt(M4*M6);
> Dxx56:=2*OM*e56*(gxx56*f1+f2)/sqrt(M5*M6);
> gxy12:=(x1-x2)*(y1-y2)/r12^2;
gyz12:=(y1-y2)*(z1-z2)/r12^2; gxy13:=(x1-x3)*(y1-y3)/r13^2; gyz13:=(y1-y3)*(z1-z3)/r13^2;
gxy14:=(x1-x4)*(y1-y4)/r14^2; gyz14:=(y1-y4)*(z1-z4)/r14^2;
gxy15:=(x1-x5)*(y1-y5)/r15^2; gyz15:=(y1-y5)*(z1-z5)/r15^2;
gxy16:=(x1-x6)*(y1-y6)/r16^2; gyz16:=(y1-y6)*(z1-z6)/r16^2;
gxz12:=(x1-x2)*(z1-z2)/r12^2; gxz13:=(x1-x3)*(z1-z3)/r13^2;
gxz14:=(x1-x4)*(z1-z4)/r14^2; gxz15:=(x1-x5)*(z1-z5)/r15^2;
gxz16:=(x1-x6)*(z1-z6)/r16^2;
> gxy23:=(x2-x3)*(y2-y3)/r23^2;
gyz23:=(y2-y3)*(z2-z3)/r23^2; gxy24:=(x2-x4)*(y2-y4)/r24^2;
gyz24:=(y2-y4)*(z2-z4)/r24^2; gxy25:=(x2-x5)*(y2-y5)/r25^2; gyz25:=(y2-y5)*(z2-z5)/r25^2;
gxy26:=(x2-x6)*(y2-y6)/r26^2; gyz26:=(y2-y6)*(z2-z6)/r26^2;
gxz23:=(x2-x3)*(z2-z3)/r23^2; gxz24:=(x2-x4)*(z2-z4)/r24^2;
gxz25:=(x2-x5)*(z2-z5)/r25^2; gxz26:=(x2-x6)*(z2-z6)/r26^2;
> gxy34:=(x3-x4)*(y3-y4)/r34^2; gyz34:=(y3-y4)*(z3-z4)/r34^2;
gxy35:=(x3-x5)*(y3-y5)/r35^2; gyz35:=(y3-y5)*(z3-z5)/r35^2;
gxy36:=(x3-x6)*(y3-y6)/r36^2; gyz36:=(y3-y6)*(z3-z6)/r36^2;
gxz34:=(x3-x4)*(z3-z4)/r34^2; gxz35:=(x3-x5)*(z3-z5)/r35^2;
gxz36:=(x3-x6)*(z3-z6)/r36^2;
> gxy45:=(x4-x5)*(y4-y5)/r45^2;
gyz45:=(y4-y5)*(z4-z5)/r45^2; gxy46:=(x4-x6)*(y4-y6)/r46^2;
gyz46:=(y4-y6)*(z4-z6)/r46^2; gxz45:=(x4-x5)*(z4-z5)/r45^2;
gxz46:=(x4-x6)*(z4-z6)/r46^2;
> gxy56:=(x5-x6)*(y5-y6)/r56^2;
gyz56:=(y5-y6)*(z5-z6)/r56^2; gxz56:=(x5-x6)*(z5-z6)/r56^2;
> Dxy12:=2*OM*e12*(gxy12*f1)/sqrt(M1*M2);
Dxy13:=2*OM*e13*(gxy13*f1)/sqrt(M1*M3); Dxy14:=2*OM*e14*(gxy14*f1)/sqrt(M1*M4);
Dxy15:=2*OM*e15*(gxy15*f1)/sqrt(M1*M5); Dxy16:=2*OM*e16*(gxy16*f1)/sqrt(M1*M6);
> Dxy23:=2*OM*e23*(gxy23*f1)/sqrt(M2*M3);
Dxy24:=2*OM*e24*(gxy24*f1)/sqrt(M2*M4); Dxy25:=2*OM*e25*(gxy25*f1)/sqrt(M2*M5);
Dxy26:=2*OM*e26*(gxy26*f1)/sqrt(M2*M6);
> Dxy34:=2*OM*e34*(gxy34*f1)/sqrt(M3*M4);
Dxy35:=2*OM*e35*(gxy35*f1)/sqrt(M3*M5); Dxy36:=2*OM*e36*(gxy36*f1)/sqrt(M3*M6);
> Dxy45:=2*OM*e45*(gxy45*f1)/sqrt(M4*M5);
Dxy46:=2*OM*e46*(gxy46*f1)/sqrt(M4*M6);
> Dxy56:=2*OM*e56*(gxy56*f1)/sqrt(M5*M6);
> Dyz12:=2*OM*e12*(gyz12*f1)/sqrt(M1*M2);
Dyz13:=2*OM*e13*(gyz13*f1)/sqrt(M1*M3); Dyz14:=2*OM*e14*(gyz14*f1)/sqrt(M1*M4);
Dyz15:=2*OM*e15*(gyz15*f1)/sqrt(M1*M5); Dyz16:=2*OM*e16*(gyz16*f1)/sqrt(M1*M6);
> Dyz23:=2*OM*e23*(gyz23*f1)/sqrt(M2*M3);
Dyz24:=2*OM*e24*(gyz24*f1)/sqrt(M2*M4); Dyz25:=2*OM*e25*(gyz25*f1)/sqrt(M2*M5);
Dyz26:=2*OM*e26*(gyz26*f1)/sqrt(M2*M6);
> Dyz34:=2*OM*e34*(gyz34*f1)/sqrt(M3*M4); Dyz35:=2*OM*e35*(gyz35*f1)/sqrt(M3*M5);
Dyz36:=2*OM*e36*(gyz36*f1)/sqrt(M3*M6);
> Dyz45:=2*OM*e45*(gyz45*f1)/sqrt(M4*M5);
Dyz46:=2*OM*e46*(gyz46*f1)/sqrt(M4*M6);
> Dyz56:=2*OM*e56*(gyz56*f1)/sqrt(M5*M6);
> Dxz12:=2*OM*e12*(gxz12*f1)/sqrt(M1*M2);
Dxz13:=2*OM*e13*(gxz13*f1)/sqrt(M1*M3); Dxz14:=2*OM*e14*(gxz14*f1)/sqrt(M1*M4);
Dxz15:=2*OM*e15*(gxz15*f1)/sqrt(M1*M5); Dxz16:=2*OM*e16*(gxz16*f1)/sqrt(M1*M6);
> Dxz23:=2*OM*e23*(gxz23*f1)/sqrt(M2*M3); Dxz24:=2*OM*e24*(gxz24*f1)/sqrt(M2*M4);
Dxz25:=2*OM*e25*(gxz25*f1)/sqrt(M2*M5); Dxz26:=2*OM*e26*(gxz26*f1)/sqrt(M2*M6);
> Dxz34:=2*OM*e34*(gxz34*f1)/sqrt(M3*M4);
Dxz35:=2*OM*e35*(gxz35*f1)/sqrt(M3*M5); Dxz36:=2*OM*e36*(gxz36*f1)/sqrt(M3*M6);
> Dxz45:=2*OM*e45*(gxz45*f1)/sqrt(M4*M5);
Dxz46:=2*OM*e46*(gxz46*f1)/sqrt(M4*M6);
> Dxz56:=2*OM*e56*(gxz56*f1)/sqrt(M5*M6);
> Dyy12:=2*OM*e12*(gyy12*f1+f2)/sqrt(M1*M2);
Dyy13:=2*OM*e13*(gyy13*f1+f2)/sqrt(M1*M3);
Dyy14:=2*OM*e14*(gyy14*f1+f2)/sqrt(M1*M4); Dyy15:=2*OM*e15*(gyy15*f1+f2)/sqrt(M1*M5);
Dyy16:=2*OM*e16*(gyy16*f1+f2)/sqrt(M1*M6);
> Dyy23:=2*OM*e23*(gyy23*f1+f2)/sqrt(M2*M3);
Dyy24:=2*OM*e24*(gyy24*f1+f2)/sqrt(M2*M4);
Dyy25:=2*OM*e25*(gyy25*f1+f2)/sqrt(M2*M5);
Dyy26:=2*OM*e26*(gyy26*f1+f2)/sqrt(M2*M6);
> Dyy34:=2*OM*e34*(gyy34*f1+f2)/sqrt(M3*M4);
Dyy35:=2*OM*e35*(gyy35*f1+f2)/sqrt(M3*M5);
Dyy36:=2*OM*e36*(gyy36*f1+f2)/sqrt(M3*M6);
> Dyy45:=2*OM*e45*(gyy45*f1+f2)/sqrt(M4*M5);
Dyy46:=2*OM*e46*(gyy46*f1+f2)/sqrt(M4*M6);
> Dyy56:=2*OM*e56*(gyy56*f1+f2)/sqrt(M5*M6);
> Dzz12:=2*OM*e12*(gzz12*f1+f2)/sqrt(M1*M2);
Dzz13:=2*OM*e13*(gzz13*f1+f2)/sqrt(M1*M3);
Dzz14:=2*OM*e14*(gzz14*f1+f2)/sqrt(M1*M4);
Dzz15:=2*OM*e15*(gzz15*f1+f2)/sqrt(M1*M5);
Dzz16:=2*OM*e16*(gzz16*f1+f2)/sqrt(M1*M6);
> Dzz23:=2*OM*e23*(gzz23*f1+f2)/sqrt(M2*M3);
Dzz24:=2*OM*e24*(gzz24*f1+f2)/sqrt(M2*M4);
Dzz25:=2*OM*e25*(gzz25*f1+f2)/sqrt(M2*M5);
Dzz26:=2*OM*e26*(gzz26*f1+f2)/sqrt(M2*M6);
> Dzz34:=2*OM*e34*(gzz34*f1+f2)/sqrt(M3*M4);
Dzz35:=2*OM*e35*(gzz35*f1+f2)/sqrt(M3*M5);
Dzz36:=2*OM*e36*(gzz36*f1+f2)/sqrt(M3*M6);
> Dzz45:=2*OM*e45*(gzz45*f1+f2)/sqrt(M4*M5);
Dzz46:=2*OM*e46*(gzz46*f1+f2)/sqrt(M4*M6);
> Dzz56:=2*OM*e56*(gzz56*f1+f2)/sqrt(M5*M6);
> qR12:=qx*(x2-x1)+qy*(y2-y1)+qz*(z2-z1);
qR13:=qx*(x3-x1)+qy*(y3-y1)+qz*(z3-z1); qR14:=qx*(x4-x1)+qy*(y4-y1)+qz*(z4-z1);
qR15:=qx*(x5-x1)+qy*(y5-y1)+qz*(z5-z1); qR16:=qx*(x6-x1)+qy*(y6-y1)+qz*(z6-z1);
> qR23:=qx*(x3-x2)+qy*(y3-y2)+qz*(z3-z2);
qR24:=qx*(x4-x2)+qy*(y4-y2)+qz*(z4-z2); qR25:=qx*(x5-x2)+qy*(y5-y2)+qz*(z5-z2);
qR26:=qx*(x6-x2)+qy*(y6-y2)+qz*(z6-z2);
> qR34:=qx*(x4-x3)+qy*(y4-y3)+qz*(z4-z3);
qR35:=qx*(x5-x3)+qy*(y5-y3)+qz*(z5-z3); qR36:=qx*(x6-x3)+qy*(y6-y3)+qz*(z6-z3);
> qR45:=qx*(x5-x4)+qy*(y5-y4)+qz*(z5-z4);
qR46:=qx*(x6-x4)+qy*(y6-y4)+qz*(z6-z4);
> qR56:=qx*(x6-x5)+qy*(y6-y5)+qz*(z6-z5);
> with(linalg):
> Dxx11:=-Dxx12-Dxx13-Dxx14-Dxx15-Dxx16;
Dxx22:=-Dxx12-Dxx23-Dxx24-Dxx25-Dxx26; Dxx33:=-Dxx13-Dxx23-Dxx34-Dxx35-Dxx36;
Dxx44:=-Dxx14-Dxx24-Dxx34-Dxx45-Dxx46; Dxx55:=-Dxx15-Dxx25-Dxx35-Dxx45-Dxx56;
Dxx66:=-Dxx16-Dxx26-Dxx36-Dxx46-Dxx56;
> DD:=matrix([ [Dxx11,Dxx12*cos(qR12),Dxx13*cos(qR13),Dxx14*cos(qR14),Dxx15*cos(qR15),Dxx16*cos(qR16)],
[Dxx12*cos(qR12),Dxx22,Dxx23*cos(qR23),Dxx24*cos(qR24),Dxx25*cos(qR25),Dxx26*cos(qR26)],
[Dxx13*cos(qR13),Dxx23*cos(qR23),Dxx33,Dxx34*cos(qR34),Dxx35*cos(qR35),Dxx36*cos(qR36)],
[Dxx14*cos(qR14),Dxx24*cos(qR24),Dxx34*cos(qR34),Dxx44,Dxx45*cos(qR45),Dxx46*cos(qR46)],
[Dxx15*cos(qR15),Dxx25*cos(qR25),Dxx35*cos(qR35),Dxx45*cos(qR45),Dxx55,Dxx56*cos(qR56)],
[Dxx16*cos(qR16),Dxx26*cos(qR26),Dxx36*cos(qR36),Dxx46*cos(qR46),Dxx56*cos(qR56),Dxx66]]);
DD :=[-.2616243284, 0, .2072378030, .05438652540, 0, 0]
[0, -3.700273954, 2.016754116, .1478108269, .9430649428, .5926440676]
[.2072378030, 2.016754116, -22.26358616, 6.499193436, .1700002483,
13.37040056]
[.05438652540, .1478108269, 6.499193436, -13.62050449, 6.740994222,
.1781194759]
[0, .9430649428, 170000248, 6.740994222, -13.87731084, 6.023251426]
[0, 5926440676, 13.37040056, 1781194759, 6.023251426, -20.16441554]
(Курсивом отмечены результаты работы программы: матричные элементы
силовой матрицы, собственные значения силовой матрицы, и частоты колебаний
решетки для шести ветвей)
> eigenvals(DD);
.93860394, -18.23995812, -13.88758592, -4.509817993,
-.3117493339, -.2100000652*10-8
>
> o1:=sqrt(abs(-36.93860394));
o2:=sqrt(abs(-18.23995812)); o3:=sqrt(abs(-13.88758592));
o4:=sqrt(abs(-4.509817993)); o5:=sqrt(abs(-.3117493339));
o6:=sqrt(abs(-.2100000652e-810-8));
o1 := 6.077713710
o2 := 4.270826398
o3 := 3.726605147
o4 := 2.123633206
o5 := .5583451745
o6 := 2.828427125
>