{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Monte Carlo with rigid cluster moves\n",
    "\n",
    "\n",
    "In this example, a low density and low temperature simulation of Lennard-Jones particles is conducted with and without rigid cluster moves and with and without single particle translations.\n",
    "A droplet is expected to form at these conditions.\n",
    "Without cluster moves, this droplet is unable to translate and rotate freely once it is formed.\n",
    "Single particle translations are still important.\n",
    "In order to obey detailed balance with rigid cluster moves, clusters cannot be created or destroyed.\n",
    "Thus, single particle translations are still required to form a droplet.\n",
    "Without single particle translations, particle randomly dispersed will not be allowed to coalesce.\n",
    "\n",
    "The three possible combinations are listed as follows:\n",
    "- sim 1 - single particle translations only\n",
    "- sim 2 - rigid cluster moves only\n",
    "- sim 3 - both"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "# Usage: /Users/user/feasst/build/bin/fst < file.txt\n",
      "FEASST version 0.25.13\n",
      "MonteCarlo\n",
      "Configuration cubic_side_length=8 particle_type=lj:/feasst/particle/lj_new.txt\n",
      "Potential Model=LennardJones\n",
      "ThermoParams beta=4 chemical_potential=1\n",
      "Metropolis\n",
      "TrialTranslate tunable_param=2\n",
      "TrialAdd particle_type=lj\n",
      "Run until_num_particles=25\n",
      "# Initializing random number generator with seed: 1749409027\n",
      "Remove name=TrialAdd\n",
      "Tune\n",
      "Log output_file=cluster1.csv trials_per_write=1000\n",
      "CPUTime output_file=cluster1_cpu.txt trials_per_write=1000\n",
      "Movie clear_file=True output_file=cluster1.xyz trials_per_write=1000\n",
      "CheckEnergy decimal_places=6 trials_per_update=1000\n",
      "Run num_trials=1e5\n",
      " \n",
      " exit: 0\n",
      "# Usage: /Users/user/feasst/build/bin/fst < file.txt\n",
      "FEASST version 0.25.13\n",
      "MonteCarlo\n",
      "Configuration cubic_side_length=8 particle_type=lj:/feasst/particle/lj_new.txt\n",
      "Potential EnergyMap=EnergyMapAll Model=LennardJones\n",
      "ThermoParams beta=4 chemical_potential=1\n",
      "Metropolis\n",
      "TrialAdd particle_type=lj\n",
      "Run until_num_particles=25\n",
      "# Initializing random number generator with seed: 1749409028\n",
      "Remove name=TrialAdd\n",
      "NeighborCriteria energy_maximum=-0.5\n",
      "TrialRigidCluster neighbor_index=0\n",
      "AnalyzeCluster output_file=cluster2_cluster.csv trials_per_write=1000\n",
      "Tune\n",
      "Log output_file=cluster2.csv trials_per_write=1000\n",
      "CPUTime output_file=cluster2_cpu.txt trials_per_write=1000\n",
      "Movie clear_file=True output_file=cluster2.xyz trials_per_write=1000\n",
      "CheckEnergy decimal_places=6 trials_per_update=1000\n",
      "Run num_trials=1e5\n",
      " \n",
      " exit: 0\n",
      "# Usage: /Users/user/feasst/build/bin/fst < file.txt\n",
      "FEASST version 0.25.13\n",
      "MonteCarlo\n",
      "Configuration cubic_side_length=8 particle_type=lj:/feasst/particle/lj_new.txt\n",
      "Potential EnergyMap=EnergyMapAll Model=LennardJones\n",
      "ThermoParams beta=4 chemical_potential=1\n",
      "Metropolis\n",
      "TrialTranslate tunable_param=2\n",
      "TrialAdd particle_type=lj\n",
      "Run until_num_particles=25\n",
      "# Initializing random number generator with seed: 1749409028\n",
      "Remove name=TrialAdd\n",
      "NeighborCriteria energy_maximum=-0.5\n",
      "TrialRigidCluster neighbor_index=0\n",
      "AnalyzeCluster output_file=cluster3_cluster.csv trials_per_write=1000\n",
      "Tune\n",
      "Log output_file=cluster3.csv trials_per_write=1000\n",
      "CPUTime output_file=cluster3_cpu.txt trials_per_write=1000\n",
      "Movie clear_file=True output_file=cluster3.xyz trials_per_write=1000\n",
      "CheckEnergy decimal_places=6 trials_per_update=1000\n",
      "Run num_trials=1e5\n",
      " \n",
      " exit: 0\n"
     ]
    }
   ],
   "source": [
    "for sim in [1, 2, 3]: # see combinations described above\n",
    "    params={\"sim\": sim,\n",
    "            \"potential\": \"Potential Model=LennardJones\",\n",
    "            \"translate\": \"TrialTranslate tunable_param=2\",\n",
    "            \"rigid\": \"\",\n",
    "            \"trials_per\": int(1e3)}\n",
    "    if sim != 1:\n",
    "        params[\"potential\"]=\"Potential Model=LennardJones EnergyMap=EnergyMapAll\"\n",
    "    if sim == 2:\n",
    "        params[\"translate\"]=\"\"\n",
    "    if sim == 2 or sim == 3:\n",
    "        params[\"rigid\"] = \"\"\"\n",
    "NeighborCriteria energy_maximum=-0.5\n",
    "TrialRigidCluster neighbor_index=0\n",
    "AnalyzeCluster trials_per_write {trials_per} output_file cluster{sim}_cluster.csv\"\"\".format(**params)\n",
    "    script=\"\"\"\n",
    "MonteCarlo\n",
    "Configuration cubic_side_length=8 particle_type=lj:/feasst/particle/lj_new.txt\n",
    "{potential}\n",
    "ThermoParams beta=4 chemical_potential=1\n",
    "Metropolis\n",
    "{translate}\n",
    "TrialAdd particle_type=lj\n",
    "Run until_num_particles=25\n",
    "Remove name=TrialAdd\n",
    "{rigid}\n",
    "Tune\n",
    "Let [write]=trials_per_write {trials_per} output_file cluster{sim}\n",
    "Log [write].csv\n",
    "CPUTime [write]_cpu.txt\n",
    "Movie [write].xyz clear_file True\n",
    "CheckEnergy trials_per_update {trials_per} decimal_places 6\n",
    "Run num_trials 1e5\n",
    "\"\"\".format(**params)\n",
    "    \n",
    "    with open('script.txt', 'w') as file: file.write(script)\n",
    "    import subprocess\n",
    "    syscode = subprocess.call(\"../../../build/bin/fst < script.txt > script.log\", shell=True, executable='/bin/bash')\n",
    "    with open('script.log', 'r') as file: print(file.read(), '\\n', 'exit:', syscode)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0    24.298629\n",
      "Name: average, dtype: float64\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "# #import math\n",
    "import pandas as pd\n",
    "\n",
    "#         end_to_end=pd.read_csv('end_to_end.txt')\n",
    "#         self.assertAlmostEqual(end_to_end['average'][0], math.sqrt(params['num_monomers']), delta=20)\n",
    "#         rg=pd.read_csv('rg.txt')\n",
    "#         self.assertAlmostEqual(rg['average'][0], params['num_monomers']/6, delta=200)\n",
    "        \n",
    "# With rigid cluster moves only, clusters cannot coalesce or break up.\n",
    "# So the clusters in the original configuration remain.\n",
    "for sim in [3]:\n",
    "    if sim == 3:\n",
    "        cluster=pd.read_csv('cluster3_cluster.csv')\n",
    "        print(cluster['average'])\n",
    "        assert cluster['average'][0] > 8"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Take a look at the cpu files. Sim 1 is faster than sim 2 because cluster moves are more expensive. Even when clusters are small, energy map book keeping and cluster calculations can be slow compared to a simple LJ potential.\n",
    "\n",
    "Sim 3 is even slower than sim 2 because sim 3 forms larger clusters. Each step moves many particles and is more expensive.\n",
    "\n",
    "Note that, even if sim 3 takes the longest time to complete a given number of steps, it can still be considered the most efficient sim due to improved sampling. Essentially, sims 1 and 2 are non ergodic and get trapped in a particular region of phase space, while sim 3 can sample all of phase space."
   ]
  }
 ],
 "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.13.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
