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

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

Ферромагнетиками называются вещества, которые могут обладать ненулевым магнитным моментом даже в отсутствие внешнего магнитного поля. Состояние ферромагнитного вещества характеризуется объемной плотностью магнитного момента $\vec{M}$, которая называется намагниченностью. Равновесному распределению намагниченности отвечает локальный минимум энергии системы. При этом ферромагнетик разделяется на области с различными ориентациями магнитного момента — домены. Домены отделены друг от друга доменными стенками — областями, в которых намагниченность поворачивается от одного направления к другому. Если энергия не минимальна, а также если ферромагнетик помещается во внешнее магнитное поле, намагниченность меняется со временем. Поведение намагниченности описывается уравнением Ландау-Лифшица-Гильберта (ЛЛГ):
$$
\frac{\partial \vec{M}}{\partial t} = -\gamma \vec{M} \times \vec{B}_{\text{eff}} + \alpha \frac{\vec{M}}{M_0} \times \frac{\partial \vec{M}}{\partial t}.
$$ Первое слагаемое описывает прецессию намагниченности в магнитном поле, а второе — потери энергии. Здесь $\vec{B}_{\text{eff}}$ — эффективное магнитное поле, действующее на магнитный момент среды в заданной точке. Оно зависит как от намагниченности в данной точке, так и от намагниченности соседних участков образца. Здесь $\gamma$ — гиромагнитное отношение, то есть отношение магнитного момента электрона к его моменту импульса, $\alpha$ — безразмерная постоянная затухания Гильберта, $M_0$ — модуль намагниченности в данной точке пространства.

Предположим, что нам известна зависимость плотности энергии (т.е. энергии единицы объема) вещества от его намагниченности $w(\vec{M})$. Тогда эффективное магнитное поле можно найти как производную \[
\vec{B}_{\text{eff}} = - \frac{\partial w}{\partial \vec{M}},
\]
где под производной по вектору имеется в виду вектор, составленный из производныx по отдельным компонентам намагниченности
\[
B_{\text{eff}, x} =- \frac{\partial w}{\partial M_x}, \quad B_{\text{eff}, y} =- \frac{\partial w}{\partial M_y}, \quad B_{\text{eff}, z} = -\frac{\partial w}{\partial M_z}.
\]
Например, если поместить вещество во внешнее магнитное поле $\vec{B}_0$, соответствующая добавка к плотности энергии будет иметь вид $$\Delta w_B = - \vec{M} \cdot \vec{B}_0,$$и соответствующее изменение эффективного магнитного поля $\Delta \vec{B}_{\text{eff}} =\vec{B}_0$.

В этой задаче вам могут потребоваться формулы для векторного произведения в декартовых координатах
\begin{align*}(\vec{a}\times\vec{b})_x& = a_y b_z -a_z b_y,\\
(\vec{a}\times\vec{b})_y& = a_z b_x -a_x b_z,\\
(\vec{a}\times\vec{b})_z& = a_x b_y -a_y b_x.
\end{align*}
Единичные вектора вдоль координатных осей $x$, $y$ и $z$ обозначим как $\vec{e}_x$, $\vec{e}_y$, $\vec{e}_z$ соответственно.

Части A, B и C задачи можно решать независимо.

Часть А. Отдельный магнитный момент (2.0 балла)

A1  ^0.30 Покажите, что если намагниченность вещества $\vec{M}$ удовлетворяет уравнению ЛЛГ, длина вектора $\vec{M}$ не меняется со временем.

Пусть намагниченность образца однородная (то есть не зависит от координат). Введем единичный вектор $\vec{m} = \vec{M}/M_0$, где $M_0 = |\vec{M}|$. Тогда уравнение ЛЛГ можно записать в виде
$$
\frac{d \vec{m}}{d t} = - \gamma \vec{m} \times \vec{B}_{\text{eff}} + \alpha \vec{m}\times\frac{d \vec{m}}{d t}.
$$ Поскольку намагниченность зависит только от времени, частную производную можно заменить на полную. Для вычислений часто удобно выразить явно производную магнитного момента по времени, в результате получится следующее уравнение
$$
\frac{d\vec{m}}{dt} = - \frac{\gamma}{1+ \alpha^2}\vec{m} \times \vec{B}_{\text{eff}} - \frac{\alpha \gamma}{1 + \alpha^2} \left((\vec{m}\cdot\vec{B}_{\text{eff}})\, \vec{m}-\vec{B}_{\text{eff}}\right).
$$

A2  ^0.90 Пусть ферромагнетик помещен во внешнее постоянное однородное магнитное поле $B_0$, направленное вдоль оси $z$. Тогда $\vec{B}_{\text{eff}} = B_0\vec{e}_z$. Найдите зависимость проекции $m_z(t)$, если в начальный момент времени $m_y = m_z = 0$, $m_x = 1$. Вам может потребоваться следующий интеграл:
$$
\int \frac{dx}{1- x^2} = \frac{1}{2} \ln \left| \frac{1+x}{1-x}\right| + C.
$$

A3  ^0.80 Покажите, что в условиях предыдущего пункта проекция магнитного момента $\vec{m}$ на плоскость $xy$, перпендикулярную магнитному полю, вращается вокруг оси $z$ c постоянной угловой скоростью. Найдите эту угловую скорость $\Omega$, выразите ответ через $\gamma$, $\alpha$, $B_0$.

Таким образом, магнитный момент прецессирует вокруг направления эффективного поля, а затухание приводит к тому, что его направление постепенно сближается с $\vec{B}_{\text{eff}}$.

Часть B. Неподвижная доменная стенка (3.5 балла)

В этой части будем рассматривать тонкую ферромагнитную проволоку с квадратным сечением, сторона квадрата равна $a$, длина проволоки $L$. Разобьем всю проволоку на $N = L/a$ кубов со стороной $a$. Будем считать, что внутри каждого куба намагниченность с достаточной точностью однородна и равна $\vec{M}_i$. Поскольку в нашей модели модуль намагниченности $M_0$ постоянен, достаточно задать единичный вектор $\vec{m}_i = \vec{M}_i/M_0$. Направим ось $x$ вдоль проволоки, а оси $y$, $z$ — перпендикулярно ей вдоль сторон сечения.

Энергия каждого кубика со стороной $a$ представляет собой сумму следующих слагаемых:

1. Энергия анизотропии, зависящая от направления магнитного момента. Для рассматриваемой проволоки $$ W_{an} = - \frac{K a^3 M_z^2}{2 M_0^2}+  \frac{K_1 a^3 M_y^2}{2 M_0^2},$$ где $K, K_1  > 0 $ — постоянные. Эта энергия минимальна, если магнитный момент направлен вдоль оси $z$.
   
2. Обменная энергия — энергия, зависящая от взаимной ориентации магнитных моментов соседних кубиков. Для одной пары она имеет вид $$W_{ex,i} = - \frac{A}{M_0^2}{\vec{M}_i} \vec{M}_{i+1}a ,$$ где постоянная $A > 0$, так что магнитным моментам энергетически выгодно ориентироваться в одну сторону. Такие вклады есть для всех пар соседних кубов.
3. Энергия взаимодействия с внешним магнитным полем $$W_B = - a^3 \vec{B}\vec{M} .$$

Обратите внимание, что здесь приведены выражения для энергии, а не плотности энергии, поэтому при вычислении компонент эффективного магнитного поля $\vec{B}_{\text{eff}}$ нужно продифференцировать энергию по соответствующей компоненте намагниченности и
разделить на
$a^3$.

Будем рассматривать вещество со следующими параметрами:

1.  $ A = 4.5 \times 10^{-10}~  Дж/м;$
2. $K = 5.0 \cdot 10^5~Дж/м^3;$
3. $\gamma  = 1.76 \cdot 10^{11} ~с ^{-1} \cdot Тл^{-1};$
4. $M_0 = 6.0 \cdot 10^5~А/м. $

Длина стороны куба равна $a = 3.0~нм$.  

B1  ^0.50 Эффективное поле, отвечающее приведенной выше энергии, можно записать в виде $$
\vec{B}_{\text{eff}, i} = B_{\text{an}} m_{i, z} \vec{e}_z + B_{\text{ex}} (\vec{m}_{i-1}+ \vec{m}_{i+1}) - k_1 B_{\text{an}}m_{i,y} \vec{e}_y + \vec{B_0}.
$$Здесь введен безразмерный численный коэффициент, задающий отношение коэффициентов анизотропии по двум осям $k_1 = K_1/K$. Везде далее считайте $k_1 = 1$. Определите формулы и численные значения для полей $B_{\text{ex}}$ и $B_{\text{an}}$. Найдите также характерную частоту прецессии магнитного момента $\omega = \gamma B_{\text{an}}$.

Далее все магнитные поля будем измерять в единицах $B_{\text{an}}$, а время — в единицах $1/\omega$, $\tau = \omega t$. Тогда обезразмеренное поле будет иметь вид ($\beta = B_{\text{ex}}/B_{\text{an}}$, $\vec{b} = \vec{B}_0/B_{\text{an}}$)
$$
\vec{b}_{\text{eff}, i} = m_{i, z} \vec{e}_z + \beta (\vec{m}_{i-1}+ \vec{m}_{i+1}) - k_1 m_{i,y} \vec{e}_y + \vec{b},
$$а уравнение ЛЛГ примет вид
$$
\frac{d \vec{m}_i}{d \tau} = - \vec{m}_i \times \vec{b}_{\text{eff}} + \alpha \vec{m}_i\times\frac{d \vec{m}_i}{d \tau}.
$$

Минимуму энергии анизотропии отвечают магнитные моменты, ориентированные параллельно оси $z$, $m_{i,z} = \pm1$. Доменная стенка между областями с $m_{i,z} = +1$ и $m_{i,z} = -1$ представляет собой переходную область, в которой магнитный момент поворачивается от положительного направления оси $z$ к отрицательному.

С помощью программы Ferromagnetic Simulator вам предстоит изучить распределение магнитных моментов в доменной стенке в зависимости от координаты $x$.

Пусть в начальный момент времени в левой половине проволоки $m_{i,z} = +1$, а в правой $m_{i,z} = -1$. Поскольку резкое изменение направления магнитного момента приводит к большой обменной энергии, эта конфигурация не оптимальна. Зададим достаточно большой коэффициент затухания $\alpha$ и рассмотрим эволюцию системы в течение некоторого времени. Она достаточно быстро окажется в состоянии с энергией, близкой к минимальной, после чего практически перестанет меняться.

В пунктах B2-B4 внешнее магнитное поле равно нулю,
$b = 0$!

B2  ^1.00 Задайте в программе значения $\beta = 100$, $b = 0$. Подберите значения времени эволюции $T$ и постоянной затухания так, чтобы к концу эволюции значения намагниченности менялись незначительно. Укажите использованные $\alpha$, $T$ в листе ответов. В папку report сохраните тот результат работы программы, который соответствует указанным $\alpha$ и $T$, назовите файл ''B2». Шириной доменной стенки $d$ будем называть расстояние между точками, где $m_z = \pm 0.5$. Получите значение $d$ в единицах $a$.

B3  ^1.50 Исследуйте зависимость ширины доменной стенки от силы обменного взаимодействия (пропорциональной величине $\beta$) в диапазоне от $\beta = 20$ до $\beta = 400$. Получите не менее 5 точек. Занесите использованные значения $\beta$ и соответствующие $d$ в Excel файле ''B3''. Не забудьте указать названия заполненных колонок. В папку report сохраните соответствующие результаты программы в файлах с названиями вида ''B3_j'' где $j$ — номер измерения. Приведите необходимые вычисления и пояснения в листах решений.

B4  ^0.50 Считая, что $d(\beta) \sim \beta^ c$, определите значение $c$ , используя результаты предыдущего пункта. Если необходимо построить график, сохраните его в Excel файле "B3".

Часть С. Движение стенки во внешнем магнитном поле (4.5 баллов)

Включим теперь внешнее магнитное поле, направленное вдоль оси $z$, $b_{0, z} = b > 0$. Доменная стенка придет в движение. Если $b < b_{\text{cr}}$, установится равномерное движение, а если $b> b_{\text{cr}}$, движение будет неравномерным. Критическое значение магнитного поля $b_{\text{cr}}$ называется пределом Уокера.
При этом распределение магнитных моментов в движущейся доменной стенке можно связать с распределением магнитных моментов в неподвижной стенке следующим образом. Пусть в отсутствии магнитного поля $m_x(x) = f(x-x_0)$, $m_y = 0$, где $x_0$ — координата середины доменной стенки. Тогда после включения магнитного поля $m_x(x) = f(x- x_0) \cos \varphi$, $m_y = f(x-x_0) \sin \varphi $, где угол $\varphi$ не зависит от координаты, а $x_0$ — координата середины доменной стенки в текущий момент времени. При равномерном движении доменной стенки угол $\varphi$ не зависит от времени.

В части C все длины измеряются в единицах $a$, все времена — в единицах
$1/\omega$, всегда
используйте значение
$\beta = 100$, $\alpha = 0.2$.

C1  ^0.20 В какую сторону будет смещаться стенка (вправо или влево)?

C2  ^2.00 Задайте в программе значения $\alpha = 0.2$, $\beta = 100$. Снимите зависимости установившейся скорости доменной стенки $v$ и угла $\varphi$ от величины внешнего поля $b$. Промеряйте диапазон от $b = 0$ до $b = b_{\text{cr}}$. Занесите использованные значения $b$ и соответствующие $v$ и $\varphi$ в Excel файл ''C2''. Не забудьте указать названия заполненных колонок. В папку report сохраните соответствующие результаты работы программы в файлах с названиями вида ''C2_j'' где $j$ — номер измерения. Считайте известным, что $b_{\text{cr}}< 0.15$. На участке $b < b_{\text{cr}}$ у вас должно получиться не менее 5 точек. Обратите внимание, что если в процессе движения центр доменной стенки окажется на расстоянии меньше 40 до одного из краев системы, решение будет существенно искажено и из него нельзя будет получить корректных результатов. Поэтому при значениях поля, близких к критическому, вам может потребоваться максимально доступная длина системы $L = 200$. При меньших значениях поля вы можете использовать меньшую длину, чтобы ускорить вычисления.

С3  ^0.40 Постройте график зависимости $v(b)$, найдите коэффициент пропорциональности $\cfrac{dv}{db}$ на линейном участке. Для построения графика используйте Excel файл ’’C2’’.

C4  ^0.60 Из теоретических рассчетов следует, что при $b<b_{\text{cr}}$ справедливо выражение $$v=v_c\sin 2\varphi.$$

Используя данные из ”C2”, постройте линеаризованный график и найдите значения $v_{\text{cr}}$ и $b_{\text{cr}}$. Для построения графика используйте Excel файл ”C2”.

C5  ^0.30 Используя численные значения из части B, выразите $v_{c}$ в м/с, а $B_{\text{cr}} = B_{\text{an}}b_{\text{cr}}$ в Тл.

C6  ^1.00 Исследуем движение при $b> b_{\text{cr}}$. При $b = 0.2$ получите зависимость координаты доменной стенки от времени. Значения координаты стенки и времени занесите в Excel файл ''C6'', не забудьте указать названия заполненных колонок. Постройте график полученной зависимости. Сохраните использованный результат работы программы в папку report, назовите файл ''C6''.

Инструкция по работе с программой Ferromagnetic Simulator

Программа Ferromagnetic Simulator численно решает уравнение ЛЛГ для системы магнитных моментов, описанной в части B. С помощью переключателя сверху вы можете выбрать часть задачи (B или С). От части задачи зависят число кубов $N$ и начальные значения магнитного момента.

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

1. Коэффициент затухания Гильберта $\alpha$ (ячейка alpha
   ), должен лежать в диапазоне $\alpha \in [0, 1]$.
   
2. Параметр $\beta$ (ячейка beta
   ), должен лежать в диапазоне $\beta \in [0, 500]$.
   
3. Внешнее магнитное поле $b$ (ячейка $b$), должно лежать в диапазоне $b \in [0, 1]$.
   
4. Время эволюции $T$ (ячейка $T$, в обезразмеренных единицах, приведенных в части B), в течение которого нужно проводить симуляцию. Тогда программа найдет решение в интервале времен $t \in [0, T]$. Время эволюции системы должно лежать в диапазоне $T \in [0, 100]$.
   

В части С вы можете дополнительно указать значение $L$ — расстояние от начального положения центра доменной стенки до правого края проволоки в единицах стороны кубика $a$. Это значение должно лежать в диапазоне $L \in [50, 200]$. Увеличение значения $L$ увеличивает время вычислений, но позволяет рассмотреть большее перемещение доменной стенки. При этом стенка всегда находится на расстоянии $L_1 = 50$ от левого края, так что общая длина системы равна $L_1  + L$.

Укажите значения параметров, после этого нажмите кнопку Solve
, и программа начнет решать уравнение. Для некоторых начальных условий это может занять несколько минут. Когда решение завершится, активизируется кнопка построения графика.

Укажите момент времени $t \in [0, T]$ и нажмите кнопку Plot. Программа
построит для момента времени $t$ сглаженные графики зависимостей $m_x(x)$, $m_y(x)$, $m_z(x)$ в одном окне. При этом фактически вычислены только значения в целочисленных точках $x = i$. Программа сохраняет значения только в 200 моментах времени, кратных $T/200$, если вы укажете другое значение, оно будет округлено до ближайшего сохраненного
.

После того, как решение получено, сохраните его с помощью кнопки Save в папку report. Пожалуйста, сохраняйте все решения, которые вы используете.
Уже сохраненное решение можно загрузить с помощью кнопки Load. После загрузки вы можете построить графики $m_x(x)$, $m_y(x)$, $m_z(x)$ для любого из 200 моментов времени, таким образом вам не требуется два раза производить вычисления для одного набора параметров.

В части B используется $N = 201$ кубик и начальные условия

1. $m_{i, z} = +1$ при $i = 0, \dots, \, 97$,
2. $m_{i, x} = 1$ при $i = 98, \dots, \, 102$,
   
3. $m_{i, z} = -1$ при $i =103, \dots, \, 200 $.
   

Остальные компоненты равны нулю.

В части C в качестве начальных условий используется вычисленное заранее статическое решение, вычисленное для доменной стенки при $\beta = 100$ для $N = 401$ кубиков. При этом лишние значения отбрасываются, чтобы соответствовать требуемой длине, заданной выше.