{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Monte Carlo Relaxation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from __future__ import print_function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import random, math"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup perfect system"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from lammps import IPyLammps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L = IPyLammps()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.units(\"lj\")\n",
    "L.atom_style(\"atomic\")\n",
    "L.atom_modify(\"map array sort\", 0, 0.0)\n",
    "\n",
    "L.dimension(2)\n",
    "\n",
    "L.lattice(\"hex\", 1.0)\n",
    "L.region(\"box block\", 0, 10, 0, 5, -0.5, 0.5)\n",
    "\n",
    "L.create_box(1, \"box\")\n",
    "L.create_atoms(1, \"box\")\n",
    "L.mass(1, 1.0)\n",
    "\n",
    "L.pair_style(\"lj/cut\", 2.5)\n",
    "L.pair_coeff(1, 1, 1.0, 1.0, 2.5)\n",
    "L.pair_modify(\"shift\", \"yes\")\n",
    "\n",
    "L.neighbor(0.3, \"bin\")\n",
    "L.neigh_modify(\"delay\", 0, \"every\", 1, \"check\", \"yes\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.image(zoom=1.6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.run(0);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "emin = L.eval(\"pe\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.dump(\"3 all movie 25 movie.mp4 type type zoom 1.6 adiam 1.0\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Disorder system"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "random.seed(27848)\n",
    "deltaperturb = 0.2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(L.system.natoms):\n",
    "    x, y = L.atoms[i].position\n",
    "    dx = deltaperturb * random.uniform(-1, 1)\n",
    "    dy = deltaperturb * random.uniform(-1, 1)\n",
    "    L.atoms[i].position = (x+dx, y+dy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.run(0);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.image(zoom=1.6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Minimize using Monte Carlo moves"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "estart = L.eval(\"pe\")\n",
    "elast = estart"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "naccept = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "energies = [estart]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "niterations = 3000\n",
    "deltamove = 0.1\n",
    "kT = 0.05"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "natoms = L.system.natoms\n",
    "\n",
    "for i in range(niterations):\n",
    "    iatom = random.randrange(0, natoms)\n",
    "    current_atom = L.atoms[iatom]\n",
    "    \n",
    "    x0, y0 = current_atom.position\n",
    "    \n",
    "    dx = deltamove * random.uniform(-1, 1)\n",
    "    dy = deltamove * random.uniform(-1, 1)\n",
    "    \n",
    "    current_atom.position = (x0+dx, y0+dy)\n",
    "    \n",
    "    L.run(1, \"pre no post no\")\n",
    "    \n",
    "    e = L.eval(\"pe\")\n",
    "    energies.append(e)\n",
    "    \n",
    "    if e <= elast:\n",
    "        naccept += 1\n",
    "        elast = e\n",
    "    elif random.random() <= math.exp(natoms*(elast-e)/kT):\n",
    "        naccept += 1\n",
    "        elast = e\n",
    "    else:\n",
    "        current_atom.position = (x0, y0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.xlabel('iteration')\n",
    "plt.ylabel('potential energy')\n",
    "plt.plot(energies)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.eval(\"pe\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "emin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "estart"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "naccept"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.image(zoom=1.6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# close dump file to access it\n",
    "L.undump(3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L.video(\"movie.mp4\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.9.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}