Monte-Carlo.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before you turn this problem in, make sure everything runs as expected. First, **restart the kernel** (in the menubar, select Kernel$\\rightarrow$Restart) and then **run all cells** (in the menubar, select Cell$\\rightarrow$Run All).\n",
    "\n",
    "Make sure you fill in any place that says `YOUR CODE HERE` or \"YOUR ANSWER HERE\", as well as your github nickname below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "NICKNAME = \"\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "eea0b4c2ba05640528fa9e4319357e1d",
     "grade": false,
     "grade_id": "cell-56cbcb736139c03b",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "# Метод Монте-Карло в физике"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "a339c3a84b55463f59280bed0c8ee43f",
     "grade": false,
     "grade_id": "cell-6576dbe1b9e30af9",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Методами Монте-Карло называют множество различных методов в совершенно разных областях от физики до экономики. Основная суть этих методов заключается в использовании случайных величин для получения различных состояний исследуемой системы. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "d7dd210ab3ea44c91133cadd14f1fa40",
     "grade": false,
     "grade_id": "cell-58e35dd9137424ae",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Одно из частых применений метода Монте-Карло в физике состоит в расчете интегралов. При этом действительно полезным и незаменимым метод становится при работе с большим количеством измерений, когда построение обычных квадратур приводит к чересчур  огромным сеткам. Например, уже в случае 10 измерений, задание обычной равномерной решетки шириной в $2^4 = 16$ узлов, приводит к общей сетке размера $2^{40} \\simeq 10^{12}$. Массив данных float такого размера будет занимать 4ТБ в памяти компьютера. Оперировать такими интегралами стандартными методами становится принципиально невозможно."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "ba6a7dc9d4cfc57bb057999fe06868da",
     "grade": false,
     "grade_id": "cell-81813c7b5296fe8f",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import scipy.special as sp\n",
    "plt.style.use('ggplot')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "408e7bacae9c68332c2d98cd2cee25eb",
     "grade": false,
     "grade_id": "cell-14d4c849aa697dc6",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "## Задание 1\n",
    "Вычислим объем многомерной сферы в большом количестве измерений."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "739ea836c17dbf030bf37ac579310071",
     "grade": false,
     "grade_id": "cell-e12d5d5a71eaa14b",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "### Задание 1.1 \n",
    "Напишите функцию, вычисляющую точный объем $d$-мерной сферы с помощью явного выражения. \n",
    "Объем многомерной сферы радиуса $R$ в $d$-измерениях дается формулой: \n",
    "$$\n",
    "V_d = \\frac{\\pi^{d/2}} {\\Gamma(1 + d/2)}\n",
    "$$\n",
    "где $\\Gamma(x)$ - гамма функция, значение которой можно получить из scipy.special."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "ac702277e841fd72344f4905bfbd2927",
     "grade": false,
     "grade_id": "cell-46a7bd98d2088625",
     "locked": false,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getSphereVolumeExact(D, R = 1.0):\n",
    "    \"\"\"\n",
    "    Функция вычисляет значение объема D-мерной сферы радиуса R по явной формуле\n",
    "\n",
    "    --------\n",
    "    Аргументы:\n",
    "    D - int, количество измерений \n",
    "    R = 1 - float, радиус сферы\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    V - float, объем сферы\n",
    "    \"\"\"\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return V"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "26cc5a1242ae5f816727e7a8a6576672",
     "grade": false,
     "grade_id": "cell-6e7a66fafbabb80c",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "1. Проверьте результат для нескольких известных вам значений. \n",
    "2. Построите график V(D)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "d54ad6aa4fa4bd1195cdd7b4080f3ff7",
     "grade": true,
     "grade_id": "cell-e3abbe900c927a0c",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "e1f704ef2557c69c68f0ba8306f0a590",
     "grade": false,
     "grade_id": "cell-6746c4683191a10a",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "### Задание 1.2\n",
    "Напишите функцию, вычисляющую объем сферы методом Монте-Карло на псевдослучайных числах (используя, например `numpy.random.rand`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "7a684a4679d07ffcedf9756be7dc7a70",
     "grade": false,
     "grade_id": "cell-1a0370e3818326f5",
     "locked": false,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getSphereVolumePseudorandom(N, D, R = 1):\n",
    "    \"\"\"\n",
    "    Функция вычисляет значение объема D-мерной сферы радиуса R методом Монте-Кало на псевдослучайных числах\n",
    "\n",
    "    --------\n",
    "    Аргументы:\n",
    "    N - int, количество случайных точек\n",
    "    D - int, количество измерений \n",
    "    R = 1 - float, радиус сферы\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    V - float, объем сферы\n",
    "    \"\"\"\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return V"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "af83d5e086486ba1208528434be0323b",
     "grade": false,
     "grade_id": "cell-470c46d55a104b03",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Постройте график зависимости относительной точности $\\epsilon$ которая достигается при фиксированном $N$ от размерности сферы $D$. Для относительной точности используйте выражение $\\epsilon = (V_{exact} - V_{random})/V_{exact}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "b316d79151ab61c96b0da29d16e959c0",
     "grade": true,
     "grade_id": "cell-bc252a77e4f5fbd7",
     "locked": false,
     "points": 2,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "64d0594d7fc0b8e470108eb2974397c3",
     "grade": false,
     "grade_id": "cell-853f758127514e1d",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "### Задание 1.3\n",
    "Напишите функцию, вычисляющую объем сферы методом Монте-Карло на квазислучайной последовательности Соболя.\n",
    "\n",
    "Для установки библиотеки `sobol_seq`:\n",
    "\n",
    "\n",
    "выполните `pip install sobol_seq` в Jupyter или\n",
    "\n",
    "выполните команду `pip install sobol_seq --user` в командном окне (если вы используете Jupyter Notebook в составе Anaconda откройте Anaconda Prompt и выполните команду pip в ней, после чего перезапустите Jupyter Notebook).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-11-17T13:15:33.843734Z",
     "start_time": "2022-11-17T13:15:33.527008Z"
    },
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "ca2bd60ff215b0c47e8fe9a1f88f2c27",
     "grade": true,
     "grade_id": "cell-cdee67dec70bf876",
     "locked": false,
     "points": 0,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "#!pip install sobol_seq <-- uncomment to install sobol_seq\n",
    "import sobol_seq\n",
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "a0708eb7467e9e14de69c825b6267544",
     "grade": false,
     "grade_id": "cell-a4e83154a99b336c",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Используйте функцию `i4_sobol_generate`. Документацию по ней можно получить выполнив следующую строчку кода:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "4dd80d4e49a11c30ff18fa6a94acc980",
     "grade": false,
     "grade_id": "cell-c2e083778576396c",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "sobol_seq.i4_sobol_generate?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "fd6332fd33d93a423b00d29d4c5ec2cf",
     "grade": false,
     "grade_id": "cell-a47b8ac1df1c3aa2",
     "locked": false,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getSphereVolumeQuasirandom(N, D, R = 1):\n",
    "    \"\"\"\n",
    "    Функция вычисляет значение объема D-мерной сферы радиуса R методом Монте-Кало на квазислучайных числах\n",
    "\n",
    "    --------\n",
    "    Аргументы:\n",
    "    N - int, количество случайных точек\n",
    "    D - int, количество измерений \n",
    "    R = 1 - float, радиус сферы\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    V - float, объем сферы\n",
    "    \"\"\"\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return V"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "e832833a2abf17f3a31f0b502d0c6269",
     "grade": false,
     "grade_id": "cell-c78358c73e3552e1",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Сравните точность вычисления объема для квазислучайной и псевдослучайной последовательностей, рассчитав  $\\epsilon(D)$, при таком же значении $N$, как в предыдущем задании. Постройте обе зависимости $\\epsilon(D)$ на одном графике в логарифмической шкале по $\\epsilon$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "cea5e5d71d2557210c2c83a3575595a3",
     "grade": true,
     "grade_id": "cell-fa2634c16c9366de",
     "locked": false,
     "points": 3,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "e383870464ab39ce5a420a82b0d9126f",
     "grade": false,
     "grade_id": "cell-e969d95dbed841c4",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Сделайте вывды о работе псевдослучайной и квазислучайной последовательностях. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "df94c0084812a2d9671fe006d76698bd",
     "grade": true,
     "grade_id": "cell-38edaea06714934d",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "source": [
    "YOUR ANSWER HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "fe4e91df7e49b802f59fe84cad36d487",
     "grade": false,
     "grade_id": "cell-5190799a80fee19b",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "# Задание 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "cc52dbc64e89d296d8859c43503ba17c",
     "grade": false,
     "grade_id": "cell-b514f20dc52808b3",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "## Задание 2.1\n",
    "Вычислите значение интеграла \n",
    "$$\n",
    "\\int d^d\\mathbf{r} \\frac{e^{-\\mathbf{r}^T\\mathbf{r}}}{1 + \\mathbf{r}^T\\mathbf{r}}, \\;\\mathbf{r} = (r_1, r_2, ..., r_d)\n",
    "$$\n",
    "по всему пространству, используя генерацию случайной величины с гауссовой плотностью.\n",
    "Проведите результаты для $d=1,2,3,4,5$, используя разумной количество точек.\n",
    "Для $d=1,2$ сравните с ответом, полученным \\textit{любым} другим способом."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "03b3f3f41ff9c708a754cad6ef9b0e0e",
     "grade": true,
     "grade_id": "cell-1a5d7fa2b658a67f",
     "locked": false,
     "points": 3,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "ccee63c4b03c4d79f02913b7313ba608",
     "grade": false,
     "grade_id": "cell-ce2fa89765e2cb96",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "## Задание 2.2\n",
    "Реализуйте функцию, возвращающую выборку случайной величины с плотностью $\\omega(x) \\propto \\sqrt{1-x^2}, x\\in[-1,1]$, используя rejection method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-11-17T13:05:20.564428Z",
     "start_time": "2022-11-17T13:05:20.558161Z"
    },
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "6e0d4d3b87a8503dafd0bb5847b17532",
     "grade": false,
     "grade_id": "cell-7e3e16764a1a98d5",
     "locked": false,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getSample(N):\n",
    "    \"\"\"\n",
    "    Функция возвращает выборку случайных величин с распределнием ~sqrt(1-x^2) \n",
    "\n",
    "    --------\n",
    "    Аргументы:\n",
    "    N - int, количество случайных точек\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    sample - numpy ndarray of float, массив случайных величин размера N\n",
    "    \"\"\"\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return sample"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "b20ba7e8d517e79fc84955c12a4f96ff",
     "grade": false,
     "grade_id": "cell-6e7e433d4b66c555",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Проверьте правильность, сгенерировав достаточно большой sample, и построив гистограмму."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "ad5194f981453ca90fb070574e726017",
     "grade": true,
     "grade_id": "cell-bf9f681eae8300e9",
     "locked": false,
     "points": 3,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "efa058240c57aa7bb6cf5420101cb234",
     "grade": false,
     "grade_id": "cell-b74adc6677b414a1",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "# Задание 3\n",
    "\n",
    "Рассмотрим решетку $n\\times n$, в которой каждой ячейке соответствует спин $\\sigma_i$ равный $+1$ или $-1$. При этом между соседними ячейками существует энергия взаимодействия\n",
    "$$\n",
    "E_i = -J \\sigma_i \\sigma_j\n",
    "$$\n",
    "Тогда энергия каждой ячейки \n",
    "$$\n",
    "E_i = -J \\sigma_i \\sum_j \\sigma_j\n",
    "$$\n",
    "где суммирование подразумевается только по соседним с $i$-ой ячейкой.\n",
    "\n",
    "При нулевой температуре энергитически более выгодным состоянием является состояние, когда все спины одинаково ориентированы. Однако при ненулевой температуре вероятность перехода в более возбужденное состояние пропорциональна $\\sim \\exp(-\\Delta E/ k T)$.\n",
    "\n",
    "Явное численное моделирование такой системы (модели Изинга) является проблематичным, поскольку количество возможных состояний системы $2^N$, где $N$ - общее число узлов, и уже для двумерной сетки $16\\times 16$ это число становится астрономически большим. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "10f464e77c2908321438cf9bdb6c8243",
     "grade": false,
     "grade_id": "cell-b3a6b8d7d9dc847b",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "На помощь приходят методы Монте-Карло. Одна из простейших его реализаций (алгоритм Метрополиса — Гастингса):\n",
    "1. Задаем случайное начальное состояние\n",
    "2. Выбираем случайную ячейку\n",
    "3. Переворачиваем её спин и считаем изменение энергии $\\delta E$ данной ячейки.\n",
    "4. Если изменение энергии отрицательное, оставляем спин перевернутым\n",
    "5. Если изменение энергии положительное, оставляем спин перевернутым с вероятностью $\\exp(-\\delta E/kT)$\n",
    "6. Повторяем пункты 2-5 нужное количество раз\n",
    "\n",
    "Спустя достаточно большое количество повторений устанавливается равновесное состояние. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "94fbd3a6424f0888f58d81a4c255fb30",
     "grade": false,
     "grade_id": "cell-1b3a73340d1b1f7e",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Реализуем данный алгоритм: \n",
    "1. Задайте начальное случайное расраспределение в виде матрицы из $+1$ и $-1$ размерности $N\\times N$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "d90cf3bba663136d8c87ae132b56cd5b",
     "grade": true,
     "grade_id": "cell-9128b635e7b2a8ed",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getInitialState(N):\n",
    "    '''\n",
    "    Функция задает случайное начальное состояние\n",
    "    ---------\n",
    "    Аргументы:\n",
    "    N - int, линейный размер решетки \n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    state - numpy ndarray of ints, массив состояния системы размера NxN\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return state"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "79b98f806bb18803e0222b4dee925fbf",
     "grade": false,
     "grade_id": "cell-e29e48763bb6018d",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "2. Напишите функцию, которая для ячейки $(i,j)$ расчитывает изменение её энергии, если изменить её знак. При этом для простоты предполагаем, во-первых, что $J=1$ (на него отнормируем kT), и во-вторых, что элементы находящиеся на верхней и нижней (левой и правой) границах находятся рядом друг с другом (топология тора).  \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "c45e1a1ce3337837eb4af34278869fb7",
     "grade": true,
     "grade_id": "cell-efa3e23518fc6507",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getDeltaE(i, j, state):\n",
    "    '''\n",
    "    Функция расчитывает и возвращает изменение энергии ячейки (i,j) в случае ее переворота (не меняя состояния state)\n",
    "    \n",
    "    ---------\n",
    "    Аргументы:\n",
    "    i - int, адресс ячейки вдоль оси 0\n",
    "    j - int, адресс ячейки вдоль оси 1\n",
    "    state - numpy ndarray of ints, массив состояния системы размера NxN\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    dE - float, изменение энергии\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return dE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "f068251f2be978138c015fefd0a03fe9",
     "grade": false,
     "grade_id": "cell-48788fd9ea8bdab0",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "3. Напишите функцию, которая $N^2$ раз для разных случайных $(i,j)$ меняет спин в зависимости от изменения энергии, согласно пунктам 4. и 5. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "689fbde57521049f5517aa6b1c75b413",
     "grade": true,
     "grade_id": "cell-609315b4a0ff1a8c",
     "locked": false,
     "points": 3,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def makeFlip(T, state):\n",
    "    '''\n",
    "    Функция N^2 раз переворачивает спин, и принимает решение оставить его перевернутым или нет, согласно правилам 4. и 5.\n",
    "    \n",
    "    ---------\n",
    "    Аргументы:\n",
    "    T - float, положительное число, безразмерный коэфициент, характеризующий температуру, равный kT/J\n",
    "    state - numpy ndarray of ints, массив состояния системы размера NxN\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    state - numpy ndarray of ints, массив нового состояния системы размера NxN\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return state"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "e23bbff7eedbee9eb3c93da644122545",
     "grade": false,
     "grade_id": "cell-69e006fc96bdf0c2",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "4. Напишите функцию, которая будет рассчитывать значение энергии всей системы (предполагаем $J=1$)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "524e8b34114ffd4c26b8adc3d516f5f6",
     "grade": true,
     "grade_id": "cell-1b5c253c53203959",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getEnergy(state):\n",
    "    '''\n",
    "    Функция, рассчитывает значение энергии всей системы\n",
    "    ---------\n",
    "    Аргументы:\n",
    "    state - numpy ndarray of ints, массив состояния системы размера NxN\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    E - float, значение энергии системы\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return E"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "8f93f892f3291b506faeeef8196ccf19",
     "grade": false,
     "grade_id": "cell-46b35827d92faf4c",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "5. Напишите функцию, которая будет рассчитывать значение намагниченности всей системы $M=\\sum \\sigma$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "ccc3dbd159e528ce8a288f9f1e1b7eb9",
     "grade": true,
     "grade_id": "cell-22bc6c63f60c8463",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "def getMagnetization(state):    \n",
    "    '''\n",
    "    Функция, рассчитывает значение намагниченности всей системы\n",
    "    ---------\n",
    "    Аргументы:\n",
    "    state - numpy ndarray of ints, массив состояния системы размера NxN\n",
    "    --------\n",
    "    Функция возвращает:\n",
    "    M - float, значение намагниченности системы\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    raise NotImplementedError()\n",
    "    return M"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "e2bf401005e850e21db5f3b53959242c",
     "grade": false,
     "grade_id": "cell-5f91305be438d924",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "В следующей ячейке кода зададим начальные параметры: \n",
    "\n",
    "**Nt** - количество точек по температуре, \n",
    "\n",
    "**N** - линейный размер решетки, \n",
    "\n",
    "**eqStep** - количество попыток переворота спина, в течение которых система приводится к равновесному сосотоянию. Обратите внимание, что это число должно увеличиваться пропорционально размеру всей сетки,\n",
    "\n",
    "**steps** - количество попыток переворота спина, в течение которых рассчитывается значения энергии и намагниченности (нужно для усреднения расчитанных значений по ансамблю состояний, близких к равновесному)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "3d08149ea9a4bb638e744e74c2916b85",
     "grade": false,
     "grade_id": "cell-79f11b7507c49140",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "N       = 10         # размер решетки NxN\n",
    "Nt      = 10         # количество точек температуры\n",
    "eqSteps = 150        # количество раз выполнения makeFlip для установления равновесия\n",
    "steps   = 30         # количество раз выполнения makeFlip для усреднения энергии и намагниченности\n",
    "\n",
    "T = np.linspace(0.5, 5, Nt); \n",
    "E, M = np.zeros(Nt), np.zeros(Nt)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "c24c0ca28646a7b6c01a814f437d8c57",
     "grade": false,
     "grade_id": "cell-81c8f6970955f96e",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "for t in range(Nt):\n",
    "    print(\"Complete\", t/Nt*100, '%\\r', end='')\n",
    "    \n",
    "    Esum = Msum = 0\n",
    "    state = getInitialState(N)\n",
    "    \n",
    "    for i in range(eqSteps):         # установление статистического равновесия\n",
    "        makeFlip(T[t], state)          \n",
    "\n",
    "    for i in range(steps):           # суммирование по разным состояниям близким к равновеснсому\n",
    "        makeFlip(T[t], state)          \n",
    "        Esum += getEnergy(state)       \n",
    "        Msum += getMagnetization(state) \n",
    "\n",
    "    E[t] = Esum / (steps * N * N)\n",
    "    M[t] = Msum / (steps * N * N)\n",
    "\n",
    "print(\"Done              \\r\", end='')\n",
    "\n",
    "_, ax = plt.subplots(1,2,figsize=(10, 5));\n",
    "\n",
    "ax[0].scatter(T, E)\n",
    "ax[0].set_xlabel(\"Temperature\");\n",
    "ax[0].set_ylabel(\"Energy \");\n",
    "\n",
    "ax[1].scatter(T, abs(M), color = 'blue')\n",
    "ax[1].set_xlabel(\"Temperature\"); \n",
    "ax[1].set_ylabel(\"Magnetization \");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "1fb7063da5118c9fe84ad7ba058f95f5",
     "grade": false,
     "grade_id": "cell-856d7c9bc2ed3bc4",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Интерпретируйте с точки зрения физики полученные графики. Описывают ли они поведение ферромагнетиков и фазовые переходы?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "29e314b20907367c58419147f6573eb4",
     "grade": true,
     "grade_id": "cell-b0c53c58b53735de",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "source": [
    "YOUR ANSWER HERE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "editable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "83e3a6d5bf845bfae57a49e8eb3e655d",
     "grade": false,
     "grade_id": "cell-a66eaf8d2e060a5a",
     "locked": true,
     "schema_version": 3,
     "solution": false,
     "task": false
    }
   },
   "source": [
    "Задайте новое сосотояние системы с достаточно большим $N$ и визуализируйте его думя цветами (например с помощью `plt.scatter`). Примените несколько раз к этому сосотоянию функциию makeFlip с низкой и с высокой температурами. Объясните наблюдаемое"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "7fe95d7a138a398d72810b8aa6fc1aa3",
     "grade": false,
     "grade_id": "cell-03dc4ad9f46d8e42",
     "locked": false,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "N = 64\n",
    "state = getInitialState(N)\n",
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "code",
     "checksum": "fc631ae525a93d79c509eb1d68d30908",
     "grade": true,
     "grade_id": "cell-48a4ce67380bd407",
     "locked": false,
     "points": 4,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "outputs": [],
   "source": [
    "plt.figure()\n",
    "# Write your code here\n",
    "makeFlip(0.5, state)\n",
    "# YOUR CODE HERE\n",
    "raise NotImplementedError()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": false,
    "nbgrader": {
     "cell_type": "markdown",
     "checksum": "9087c7c2f501fcdc5b5fb1d4b67c05d7",
     "grade": true,
     "grade_id": "cell-72a37af61409d074",
     "locked": false,
     "points": 1,
     "schema_version": 3,
     "solution": true,
     "task": false
    }
   },
   "source": [
    "YOUR ANSWER HERE"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.13"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": true,
   "autocomplete": false,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}