Модель Изинга Текст научной статьи по специальности «Физика»

МЕТОДИЧЕСКИЕ МАТЕРИАЛЫ

Вестник Сыктывкарского университета. Серия 1: Математика. Механика. Информатика. Выпуск 1 (22). 2017

УДК 537.622

МОДЕЛЬ ИЗИНГА

В. А. Устюгов

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

Ключевые слова: ферромагнетизм, модель Изинга, термодинамика.

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

Описывая ферромагнетик, мы будем предполагать, что спины расположены в ячейках квадратной решетки. Проекция каждого спина на выделенную ось может принимать всего два значения: в = ±1, т. е. спин может быть направлен вверх или вниз. Эти спины взаимодействуют с другими спинами решетки, расположенными в соседних ячейках. Согласно модели Изинга энергия системы спинов рассчитывается следующим образом:

© Устюгов В. А., 2017.

где суммирование происходит по всем парам соседних спинов, 3 > 0 — обменный интеграл [4]. Будем считать, что конфигурация спинов системы определяет ее состояние, тогда приведенная формула позволяет вычислить энергию этого состояния. Также из формулы (1) следует, что если соседние спины параллельны друг другу, то энергия их пары равна Е = — 3 < 0, откуда следует, что такое расположение является энергетически более выгодным по сравнению с антипараллельным. Следовательно, в рассматриваемой системе в отсутствие внешнего магнитного поля может наблюдаться ненулевой макроскопический магнитный момент, т. е. спонтанная намагниченность.

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

Ра - е-Е«/квт, (2)

где Еа — энергия системы в состоянии а, к в — постоянная Больцмана. Если квадратная решетка содержит N ячеек, то число возможных состояний спиновой системы равно 2м.

1. Теория усредненного поля

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

М = £ (вг) = N (вг), (3)

г

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

Добавим к рассматриваемой модели внешнее магнитное поле с напряженностью Н. Тогда энергия системы выразится следующим образом:

Е = -^^в,-в* - (4)

<У> г

где ^ — магнитный момент одного спина. В этом случае спины будут стремиться ориентироваться вдоль направления поля, так как это приводит к уменьшению энергии. Предположим теперь, что система состоит из одного спина вг = ±1. Энергия системы в этом случае равна Е± = ^^Н. Вероятности найти систему в двух ее состояниях равны:

Р+ = Ое+^Н/кв т, (5)

р_ = Се-^Н/кв Т,

где нормировочный коэффициент с учетом требования равенства единицы суммы вероятностей Р+ и Р_ вычисляется как

С = е+^Н/квТ + е_^Н/квТ ' (6)

Итак, среднее равновесное значение (вг) в этом случае можно найти как

(вг) = X] э^Р± = Р+ + Р_ = tanh{^Н/квТ). (7)

«¿=±1

Используем полученную формулу для нахождения приближенного решения задачи для ансамбля из N спинов. Приближение усредненного поля предполагает, что взаимодействие спина вг с соседними спинами эквивалентно его движению в некотором эффективном поле Heff. Заменив в формуле (7) переменную Н на Heff, мы можем вычислить значение (вг). Оценим величину Heff.

Перепишем уравнение (4) в более удобном виде:

Е = - вг - ^Н^вг, (8)

\ Чз> ) г

откуда ясно, что множитель при вг может быть представлен как ^^в^- = цHeff. Теперь предположим, что переменные в^- в последнем выражении могут быть заменены на их равновесные средние значения (в). Учитывая, что все спины бесконечно большой системы имеют одинаковую равновесную среднюю ориентацию, и полагая внешнее поле

нулевым, получим:

Яе// = 3 £(в) = 3(в), (9)

где г — число ближайших соседей спина. Учитывая (7), получаем:

(в) = 1апЬ(г3 (в)/кв Т). (10)

Мы получили уравнение, не разрешенное относительно переменной в (трансцендентное), найти решение аналитически можно только при малом значении (в). Оценим решение уравнения качественно. Изобразив на одной координатной плоскости функции, соответствующие левой и правой части уравнения (10), можно увидеть, что при достаточно низкой температуре графики пересекаются в трех точках: (в) = -во, 0,в0. Нулевое решение соответствует состоянию системы с нулевой спонтанной намагниченностью, т. е. парамагнитному состоянию. Два других решения отвечают положениям спина вверх и вниз при ферромагнитном состоянии системы. При этом в отсутствие внешнего поля эти две ориентации спина являются равновероятными.

Для нахождения корней уравнения (10) можно воспользоваться численными методами поиска нулей функции применительно к выражению

/ ((в)) = (в)- 1апЬ(г3 (в)/кв Т). (11)

Находя эти корни, мы обнаружим, что при низкой температуре система находится в ферромагнитном состоянии и спонтанная намагниченность М > 0, поскольку она пропорциональна (в). При высокой температуре система пребывает в парамагнитном состоянии. В рамках рассматриваемой теории усредненного поля переход происходит при критической температуре Тс = 4, где температура выражена в единицах 3/кв, при этом график зависимости М (Т) вблизи точки Тс имеет крутой спад. Иными словами, производная йМ/йТ обращается в бесконечность при Т ^ Тс.

Точнее определить характер зависимости вблизи критической температуры можно, разлагая гиперболический тангенс в ряд для малого аргумента: 1апЬ(х) ~ х — х3/3. Поскольку значение (в) мало, то

(в)« ^М — 1 ( У . (12)

К ' кв Т 3 V кв Т ) У '

Уравнение (12) имеет тривиальное решение (в) = 0 и

3 (квТ\ ( г3

1/2

ТЗ (кВТ—Т> -<Тс—Т>'5• <13)

где Тс = г3/кв, а в = 1/2. Параметр в называется критическим индексом. В дальнейшем покажем, что полученное значение в не является корректным, поскольку модель усредненного поля не позволяет получить корректную форму выражения для (в).

2. Метод Монте-Карло

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

Случайным или некоторым систематическим образом выбирается спин и рассчитывается изменение энергии системы Е/ в случае, если он перевернется. Если Е/ < 0, т.е. энергия системы уменьшается, то спин переворачивается, переводя систему в новое микросостояние. Чтобы система не застревала в области локального минимума энергии, требуется разрешить переходы, увеличивающие энергию, но с некоторой небольшой вероятностью. Так, если Е/ > 0, происходит следующий выбор. Генерируется случайное число с однородным распределением на отрезке [0; 1] и сравнивается с величиной ехр(—Е//кьТ). Если она больше, чем сгенерированное случайное число, то спин переворачивается, иначе система остается без изменений. Таким образом производится один шаг Монте-Карло для одного спина.

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

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

верхней и нижней границ решетки, а также правой и левой.

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

Рис. 1. Зависимость намагниченности системы из 15 х 15 спинов от времени при температурах: а) Т = 1.0; Ь) Т = 2.0; с) Т = 2.25; а) Т = 3.0

При температуре Т = 2.25 мы можем наблюдать флуктуации намагниченности системы, составляющие порядка 80 % намагниченности насыщения. Это означает, что происходит перемагничивание всей систе-

12 3

Тетрегашге

а)

Ь)

Рис. 2. Результаты моделирования: зависимость от температуры а) намагниченности системы спинов и Ь) энергии системы, приходящейся на спин

250

200

150

100

50

0

4

0

3

4

Тетрегашге

мы и свидетельствует о предельной близости фазового перехода. Точный расчет для модели Изинга дает значение критической температуры Тс = 2.27, так что при Т = 3.0 система находится в парамагнитном состоянии.

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

3. Термодинамические параметры системы спинов

Рассмотрим зависимость средней энергии спина от температуры. При низкой температуре, когда все спины ориентированы в одном направлении, каждый из четырех соседей дает вклад — 3 в энергию спина, так что полная средняя энергия системы составит — 4Ж3/2, где деление на 2 необходимо, поскольку при расчете полной энергии энергия каждой пары соседних спинов войдет в итоговое выражение дважды. При высокой температуре мы ожидаем, что в среднем из четырех соседей два будут ориентированы в одном направлении, другие два — в другом, и средняя полная энергия системы будет стремиться к нулю. Однако моделирование показывает, что энергия (Е) остается ненулевой даже при температуре Т > Тс (рис. 2Ь). Это говорит о том, что при температуре выше критической некоторые спины расположены не случайным образом (т.е. ориентации соседних спинов коррелируют), несмотря на то, что суммарная намагниченность системы равна в этом случае нулю.

По графику зависимости (Е) можно видеть, что точка критической температуры является для этой функции точкой перегиба, причем первая производная функции обращается в бесконечность. Однако функция С = ^(Е)/^Т является удельной теплоемкостью тела, которая бесконечно растет при критической температуре [2]. Теплоемкость можно представить в виде степенной зависимости как

С ~ |Т - Тс|_а, (14)

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

(Е) = Еа, (15)

а

где Еа — энергия системы в состоянии а. Аналогично вычисляется величина

(Е2) = т^Е Е2. (16)

Nm

а

Согласно флуктуационно-диссипационной теореме удельную теплоемкость можно найти как

С ■ <17)

где (ДЕ)2 = (Е2) - (Е)2. Построив график зависимости, можно обнаружить на нем ярко выраженный пик в критической точке. Это связано, как было сказано, с повышенной чувствительностью системы к внешним воздействиям при критической температуре.

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

/ (¿) = (вовг), (18)

Рис. 3. Функция корреляции системы спинов при различных температурах

где г — расстояние между спинами, выраженное в количестве узлов решетки между в0 и вг. Оно выражается целым числом и для квадратной решетки спинов размером Ь х Ь с периодическими граничными условиями не превышает Ь/2.

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

При температуре Т = 2.0 дальность влияния спина несколько больше.

Вблизи критической температуры (Т = 2.5) степень влияния спина на ближние существенно выше, чем на дальние. При дальнейшем повышении температуры величина корреляции снова уменьшается.

Такое поведение функции корреляции позволяет объяснить бесконечные значения удельной теплоемкости и магнитной восприимчиво-

сти при критической температуре. При температуре T > Tc функция с расстоянием убывает экспоненциально:

f (i) - Ci + C2exp(-ri/£), (19)

где C1 = lim f(i), ri — расстояние между спинами £ — радиус корреляции. При r¿ ^ то спины so и si не коррелируют, т. е.

f (i) ^(so)(Si) = (s)2. (20)

Таким образом, становится ясно, что C1 = (s)2, и уравнение (19) можно переписать как

f (i) - Ci = ((so - (s))(si - (s))) - C2exp(-ri/£). (21)

Таким образом, £ позволяет оценить, на каком расстоянии корреляция флуктуаций спинов s0 и si обращается в ноль. Величина £ растет при приближении к критической температуре, и точный расчет показывает, что для £ подчиняется степенному закону

£ -|T - Tc|-V, (22)

где v — также критический индекс. При T = Tc £ = то, это означает, что спиновые флуктуации некоторого выбранного спина вызывают флуктуации спинов, которые находятся на значительном расстоянии от него. Таким образом, становится понятно, что изменение состояния каждого спина будет влиять на состояние всех других спинов системы. Эта исключительная чувствительность и объясняет резкий скачок магнитной восприимчивости и удельной теплоемкости. Отметим, что описанный эффект абсолютно не учитывается в теории усредненного поля и в этом состоит причина неверной оценки в рамках теории критической температуры.

4. Заключение

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

Список литературы

1. Giordano N. J., Nakanishi H. Computational physics. Pearson/Prentice Hall, 2006. 544 p.

2. Coey J. Magnetism and Magnetic Materials. Cambridge University Press, 2010. 633 p.

3. Биндер К., Хеерманн Д. В. Моделирование методом Монте-Карло в статистической физике. М.: ФИЗМАТЛИТ, 1995. 144 с.

4. Гулд Х., Тобочник Я. Компьютерное моделирование в физике. М.: Мир, 1990. 400 с.

СГУ им. Питирима Сорокина Поступила 01.03.2017

Summary

Ustyugov V. A. The Ising model

The article provides an overview of the mathematical method of the Ising model and average field theory. We compared the values of the critical temperature, analytically derived based on the mean field theory and by numerical simulation. The causes differences of these values are discussed. Keywords: ferromagnetism, Ising model, thermodynamics.

References

1. Giordano N. J., Nakanishi H. Computational physics, Pearson/ Prentice Hall, 2006, 544 p.

2. Coey, J. Magnetism and Magnetic Materials, Cambridge University Press, 2010, 633 p.

3. Binder K., Heermann D. W. Monte Carlo methods in Statistical Physics, M.: FIZMATLIT, 1995, 144 p.

4. Gould H., Tobochnik J. Computer modelling in physics, M.: Mir, 1990, 400 p.

Для цитирования: Устюгов В. А. Модель Изинга // Вестник Сыктывкарского университета. Сер. 1: Математика. Механика. Информатика. 2017. Вып. 1 (22). C. 61-71.

For citation: Ustyugov V. A. The Ising model, Bulletin of Syktyvkar University, Series 1: Mathematics. Mechanics. Informatics, 2017, №1 (22), pp. 61-71.