Лабораторная работа №1. Метод Монте-Карло
Методы Монте-Карло — это широкий класс вычислительных алгоритмов, которые полагаются на повторной случайной выборки для получения численных результатов. Как правило, один эксперимент повторяется много раз, чтобы получить распределение вероятности неизвестной величины. Название происходит от сходства техники и игрового процесса в рулетку в реальном казино. Данные методы часто используются для решения физических и математических задач и являются наиболее полезными, когда трудно или невозможно получить точное решение. Методы Монте-Карло в основном используются в трех различных классах задач: оптимизации, численного интегрирования и генерации распределения вероятностей. Для решения физических проблем, методы Монте-Карло являются весьма полезными для моделирования систем со многими степеней свободы, например, жидкостей, неупорядоченных материалов, сильно связанных твердых тел и решетчатых структур.
Монте-Карло моделирование с помощью алгоритма Метрополиса зарекомендовало себя как мощный инструмент для систематического изучения магнитных свойств наночастиц и их ансамблей. Два основных преимущества этой методики:
1. изучение наночастиц возможно вплоть до атомных масштабов, поэтому детали их структуры могут быть изучены;
2. корректный учет температуры при проведении численных экспериментов.
Алгоритм Метрополиса состоит в следующем:
1) Выбирается случайный спин в двумерном массиве и вычисляется энергия его взаимодействия с его соседями, т.н. энергия «старой» или исходной конфигурации: \(E_1\). После этого данный спин переворачивается и вычисляется энергия «новой» конфигурации: \(E_2\).
2) Энергию «новой» конфигурации сравниваем с энергией «старой» конфигурации. «Новая» конфигурация принимается и становится исходной для следующего шага, если \(E_1>E_2\), в противном случае рассчитывается вероятность переворота данного спина \(p(E_1\rightarrow E_2)\):
\begin{equation}
p(E_1 \rightarrow E_2)=\left\{ \begin{aligned}
1 \text{, если } E_1>E_2 \\
e^\frac{E_2-E_1}{T} \text{, если } E_2\leq E_1
\end{aligned} \right\}
\end{equation}
Нужно создать 2 цикла for: первый цикл по температуре, второй по МК шагам. Необходимо выполнить 1000000 МК шагов для каждого шага по температуре. Цикл по температуре в интервале от 0.01 до 4, с шагом 0.1
Модель Изинга
Модель Изинга является одной из самых простых математических моделей ферромагнетизма в статистической механике, допускающих точной решение. Математическая модель рассматривает совокупность дискретных переменных (значения магнитных моментов атомных спинов), которые могут принимать одно из двух значений, соответствующих одному из двух состояний. Для наглядности спины обычно располагаются в виде графа. Модель позволяет дать наглядную интерпретацию фазового перехода II рода. Двумерная квадратная решетка спинов Изинга является одной из самых простых статистических моделей, позволяющей строгое описание фазового перехода.
Рассмотрим множество смежных узлов, которые формируют двумерную решетку. Для каждого узла j соответствует дискретная переменная \(S_j\) принимающая значения +1 или −1. У каждого спина есть 4 соседа (справа, слева, сверху и снизу). Энергия их взаимодействия рассчитывается по формуле
\begin{equation*}
E=-J\sum_{j=1}^4 S_i S_j
\end{equation*}
Полная энергия системы учитывает парные взаимодействия спинов i, j, интенсивность взаимодействия определяется константой J=1, то энергия одного из \(2^N\) возможных состояний системы задается уравнением вида:
\begin{equation*}
E=-J\frac{\sum_{i\neq j=1}S_i S_j}{N},
\end{equation*}
где N — размер системы.
Под конфигурацией понимается магнитное состояние системы спинов, характеризующиеся уникальным распределением значений \(S_i\).
Если для всех пар i, j:
1. J >0 то взаимодействие называется «ферромагнитным»;
2. J <0 то взаимодействие называется «антиферромагнитным»;
3. J =0 то спины не взаимодействующие.
Задание
Язык разработки желателен С/С++, но Вы можете использовать любой на Ваше усмотрение. Все расчеты необходимо реализовывать через функции, наче Ваша программа не будет засчитана.
1. Реализовать функцию создания и заполнения двумерного массива, значениями +1 и −1 случайным образом. Размер массива 50×50.
2. В main() реализовать генерацию номера случайного спина: необходимо через rand() генерировать случайным образом индексы элемента Вашего двумерного массива спинов: i, j : например: array[i][j] = 1 или array[i][j] = -1, где значения +1 или −1 – это значения спина \(S_i\) .
3. Реализовать определение соседей у случайно выбранных спинов. А также учет граничных условий: например, если выбранный спин стоит в верхней строке и по сути не имеет реального соседа сверху, то за верхнего соседа принимается сосед из того же столбца, но с самого нижнего ряда. По аналогии для других соседей.
4. Реализовать функцию расчета энергии для одного спина – в эту функцию (при вызове ее в main()) необходимо передавать номер случайно выбранного спина \(S_i\) (его индексы i, j ) и производить расчет энергии взаимодействия этого спина с его 4 соседями, и функцию для расчета энергии всей системы.
5. Реализовать функцию расчета полной энергии системы и средней энергии на 1 спин.
6. Реализовать функцию расчета намагниченности М всей системы по формуле:
\begin{equation*}
M=\frac{\sum_{i=1}^N S_i}{N}
\end{equation*}
7. Рассчитать теплоемкость системы по формуле:
\begin{equation*}
C=\frac{<E^2>-<E>^2}{T^2}
\end{equation*}
8. Рассчитать среднеквадратичную ошибку эксперимента:
\begin{equation*}
\sqrt{\frac{1}{n-1}\sum_{i-1}^n (x_i-\overline{x})^2}
\end{equation*}
полученная величина будет характеризовать ошибку эксперимента для каждой точки на графике в виде ± ошибки (yerrorbars в gnuplot). Необходимо рассчитать ошибку для энергии, намагниченности и теплоемкости.
9. Организовать вывод на экран элементов массива, значения энергии одного случайно выбранного спина, энергии всей системы и намагниченности всей системы.
10. Используя gnuplot построить графики функций E(T), M(T), C(T) с учетом среднеквадратичной ошибки.