spreading_algorithm_and_code.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Algorithm and (Jupyter) python code"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### We start by importing the necessary libraries and widgets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "4eabba64bb044936879300a757c6b9c9",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "VBox(children=(HBox(children=(Output(), Output())), HBox(children=(Button(description='Next Spread', style=But…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "%config InlineBackend.close_figures=False;\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "import ipywidgets as widgets\n",
    "import numpy as np\n",
    "from matplotlib import rcParams;\n",
    "import ipywidgets as widgets;\n",
    "from ipywidgets import Layout;\n",
    "from IPython.display import display,clear_output;\n",
    "import pylab as lp;\n",
    "from pylab import *;\n",
    "from scipy.ndimage import measurements;\n",
    "from IPython.display import display,clear_output;\n",
    "import warnings;\n",
    "warnings.filterwarnings('ignore');\n",
    "import networkx as nx\n",
    "import random\n",
    "import os\n",
    "import scipy.io as sio\n",
    "from scipy.io import loadmat\n",
    "from pathlib import Path\n",
    "import itertools\n",
    "import collections\n",
    "import ipywidgets as widgets;\n",
    "from ipywidgets import Layout\n",
    "from IPython.display import display,clear_output\n",
    "from IPython.core.display import display, HTML\n",
    "from pathlib import Path\n",
    "from matplotlib import colors\n",
    "from IPython.display import HTML as html_print\n",
    "from IPython.display import display\n",
    "\n",
    "out1 = widgets.Output()\n",
    "out2 = widgets.Output()\n",
    "\n",
    "button=widgets.Button(description='Next Spread')\n",
    "buttonres=widgets.Button(description='Reset')\n",
    "hbox=widgets.HBox(children=(button,buttonres))\n",
    "vbox1=widgets.HBox(children=(out1,out2))\n",
    "vbox=widgets.VBox(children=(vbox1,hbox))\n",
    "display(vbox)\n",
    "matplotlib.rcParams['figure.figsize'] = [6, 6]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Using python networkx library, we create a graph g, associated with a 2D square grid (side N=20). We define $n_{ng}$ as a list of the nodes grid and $n_{inf}$ the number of infected nodes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "N=20;\n",
    "ninf=int(2.5*N);\n",
    "#############################################\n",
    "g = nx.grid_2d_graph(N, N);\n",
    "nng=list(g.nodes());"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Using python networkx library, we create a graph g, associated with a 2D square grid (side N=20). We define $n_{ng}$ as a list of the nodes grid and $n_{inf}$ the number of infected nodes. We also define $infected$, a distribution that defines the configuration of infected squares. Three distributions are suggested: 1-random infection; 2-diagonal infection; 3- two-sides infection. By default our code proposes the 3rd distribution (of infected squares). For a diagonal or random distribution, just uncommente the associated block, and disable (comment) the two-sides infection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [],
   "source": [
    "##############################################\n",
    "###three distribution-types of infections\n",
    "##############################################\n",
    "###1-random infection\n",
    "#infected=random.sample(nng,ninf)\n",
    "###end of random infection distrb.\n",
    "##############################################\n",
    "###2-diagonal infection\n",
    "#diagonal = np.arange(N-1,N*N-1,N-1);\n",
    "#distrib=diagonal;\n",
    "#infitems=[]\n",
    "#for i in range(len(distrib)):\n",
    "#    infitems.append([n for n in nng[distrib[i]]])\n",
    "#infected=list(tuple(map(tuple, infitems)))\n",
    "###end of diagonal infection distrb.\n",
    "##############################################\n",
    "###3-two-sides infection\n",
    "xside=np.arange(0,N*N-1,2*N)\n",
    "yside=np.arange(N*N-N+1,N*N,2)\n",
    "twosides=np.concatenate((xside,yside), axis=None)\n",
    "distrib=twosides;\n",
    "infitems=[]\n",
    "for i in range(len(distrib)):\n",
    "    infitems.append([n for n in nng[distrib[i]]])\n",
    "infected=list(tuple(map(tuple, infitems)))\n",
    "###end of two-sides infection distrb.\n",
    "##############################################\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ismember(A,B):\n",
    "    U, ixU = np.unique(A , return_index=True);\n",
    "    V = np.in1d(U, B, assume_unique=True);\n",
    "    ismember = V*1;\n",
    "    return ismember;\n",
    "##############################################"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### We define a spreading function $spreading(b)$ that spreads infection from the step n to the step n+1. This function creates a new graph $G_{inf}$ in which all infected nodes are stored computes $allneighbors$ function that finds all neighbors of infected squares. $neighb2times$ is a new function that counts which nodes are neighbors of infected squares more than once. We also need to compute $perim$ function that computes the perimeter of the whole infected area. At the end of this code snippet, we create $data\\_saver.mat$, a .mat file to store data associated with the perimeter of infected area. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [],
   "source": [
    "step=0;\n",
    "def spreading(b):\n",
    "    global step;    \n",
    "    step +=1; \n",
    "    M=[]\n",
    "    if Path('G_inf.gpickle').is_file():\n",
    "        G_inf = nx.read_gpickle(\"G_inf.gpickle\")\n",
    "    else:\n",
    "        G_inf=nx.Graph()\n",
    "        G_inf.add_nodes_from(infected)\n",
    "    inf_nodes=G_inf.nodes();\n",
    "    Ninf=list(inf_nodes)\n",
    "    NbInf=len(Ninf);\n",
    "    inf_array=np.asarray(Ninf)\n",
    "    \n",
    "    arrayinf=np.array([list(inf_array[:,0]),list(inf_array[:,1])])\n",
    "     \n",
    "    inf_index=np.ravel_multi_index(arrayinf, (N,N))\n",
    "    lix=np.arange(N*N)\n",
    "    izInf=ismember(lix,inf_index)\n",
    "    inf_nodes_inmat0=izInf*(1+lix)\n",
    "    inf_nodes_inmat1=inf_nodes_inmat0.reshape(N,N)\n",
    "    inf_nodes_inmat=np.where(inf_nodes_inmat1>0,1,0)\n",
    "        \n",
    "    Ninf=list(G_inf.nodes())\n",
    "    NbInf=len(Ninf);\n",
    "    \n",
    "    neighb=[]\n",
    "    for i in range(len(Ninf)):\n",
    "        neighb.append([n for n in g.neighbors(Ninf[i])])\n",
    "    \n",
    "    allneighbors0=list(itertools.chain.from_iterable(neighb));\n",
    "    allneighbors= [i for i in allneighbors0 if i not in Ninf ];\n",
    "    Nb_allneighbors= len(allneighbors)\n",
    "    neighb2times=[item for item,\n",
    "                 count in collections.Counter(allneighbors).items()\n",
    "                 if count > 1]\n",
    "    NodesContam = [i for i in neighb2times if i not in Ninf ]   \n",
    "    perim=Nb_allneighbors-len(NodesContam);   # perim function computes the perimeter of infected area                 \n",
    "\n",
    "    b=.7\n",
    "    if len(neighb2times)>0:\n",
    "        cont_array=np.asarray(NodesContam)\n",
    "        cmap = colors.ListedColormap(['white','black',[b,b,b]])\n",
    "    else:\n",
    "        cont_array=inf_array\n",
    "        cmap = colors.ListedColormap(['cyan',[b,b,b]])\n",
    "    arraycont=np.array([list(cont_array[:,0]),list(cont_array[:,1])]);\n",
    "    cont_index=np.ravel_multi_index(arraycont, (N,N));\n",
    "    lix=np.arange(N*N)\n",
    "\n",
    "    izCont=ismember(lix,cont_index)\n",
    "    cont_nodes_inmat0=izCont*(1+lix)\n",
    "    cont_nodes_inmat1=cont_nodes_inmat0.reshape(N,N)\n",
    "    cont_nodes_inmat=np.where(cont_nodes_inmat1>0,2,0)\n",
    "    \n",
    "    with out1:\n",
    "        fig1, axes1 = plt.subplots();\n",
    "        img=inf_nodes_inmat.transpose()+cont_nodes_inmat.transpose()\n",
    "        Ninf=list(inf_nodes)\n",
    "        NbInf=len(Ninf);\n",
    "        im = axes1.imshow(img, cmap);\n",
    "        axes1.set_xticks(np.arange(-.5, N, 1));\n",
    "        axes1.set_yticks(np.arange(-.5, N, 1));\n",
    "        axes1.xaxis.set_tick_params(width=0);\n",
    "        axes1.yaxis.set_tick_params(width=0);\n",
    "        axes1.set_xlim([-1,N]);\n",
    "        axes1.set_ylim([-1,N]);\n",
    "        axes1.grid(linestyle=\"-\", color='black');\n",
    "        axes1.set_yticklabels([])\n",
    "        axes1.set_xticklabels([]) \n",
    "        plt.title('pemimeter= '+str(perim)+', Nb. of inf. squ.= '+str(NbInf)+\n",
    "                  ', Nb. of cont. squ.= '+str(len(NodesContam)))\n",
    "        plt.show(fig1)\n",
    "\n",
    "        clear_output(wait=True)\n",
    "    \n",
    "    for i in range(len(neighb2times)):\n",
    "        Ninf.append(tuple([n for n in neighb2times[i]]))\n",
    "        G_inf.add_nodes_from(Ninf)\n",
    "        nx.write_gpickle(G_inf,\"G_inf.gpickle\")\n",
    "\n",
    "    Ninf=list(G_inf.nodes());\n",
    "    NbInf2=len(Ninf);\n",
    "    perim2=Nb_allneighbors-len(NodesContam);\n",
    "    Nbcont2= len(NodesContam)\n",
    "    Nbcont= len(NodesContam)\n",
    "\n",
    "    \n",
    "    if Path('data_saver.mat').is_file():\n",
    "        datamat=loadmat('data_saver.mat')\n",
    "        perim0=datamat['perim']\n",
    "        perim=perim0[0]\n",
    "        perim=np.concatenate((perim,perim2), axis=None)\n",
    "        NbInf0=datamat['NbInf']\n",
    "        NbInf=NbInf0[0]\n",
    "        NbInf=np.concatenate((NbInf,NbInf2), axis=None)\n",
    "        Nbcont0=datamat['Nbcont']\n",
    "        Nbcont=Nbcont0[0]\n",
    "        Nbcont=np.concatenate((Nbcont,Nbcont2), axis=None)\n",
    "        sio.savemat('data_saver.mat', {'NbInf':NbInf,'Nbcont':Nbcont,'perim':perim})\n",
    "    else:\n",
    "        sio.savemat('data_saver.mat', {'NbInf':NbInf,'Nbcont':Nbcont,'perim':perim})\n",
    "        \n",
    "\n",
    "    datamat=loadmat('data_saver.mat')\n",
    "    perim0=datamat['perim']\n",
    "    perim=perim0[0]\n",
    "    NbInf0=datamat['NbInf']\n",
    "    NbInf=NbInf0[0]\n",
    "    xx=np.arange(len(perim))\n",
    " \n",
    "    with out2:\n",
    "        fig2, axes2 = plt.subplots()\n",
    "        clear_output(wait=True)\n",
    "        axes2.plot(xx,NbInf/(N*N),'ro-');\n",
    "        axes2.plot(xx,10*Nbcont/(N*N),'go-');\n",
    "        axes2.plot(xx,(N*N-NbInf)/(N*N),'bo-');#(perim[0])\n",
    "        axes2.set_xticks(np.arange(-0.2, 3*N, N/2));\n",
    "        axes2.set_yticks(np.arange(-0.2, 3, .2));\n",
    "        axes2.set_xlim([-.2,3*N]);\n",
    "        axes2.set_ylim([-.2,3]);\n",
    "        axes2.grid();\n",
    "        plt.title('Spreading functions: healthy(b)/infected (r)/contamin.(g)');\n",
    "        axes2.set_xlabel('Contamination step');\n",
    "        axes2.set_ylabel('NbInf, perim');\n",
    "\n",
    "        plt.show(fig2)\n",
    "        clear_output(wait=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### We end our program by setting up the functions of our widgets: new spreading configuration and a reset operation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "button.on_click(spreading);\n",
    "spreading(None);\n",
    "\n",
    "def cstr(s, color='black'):\n",
    "    return \"<text style=color:{}>{}</text>\".format(color, s)\n",
    "\n",
    "def print_color(t):\n",
    "    display(html_print(' '.join([cstr(ti, color=ci) for ti,ci in t])))\n",
    "    \n",
    "\n",
    "def reset(b):\n",
    "    if os.path.exists(\"G_inf.gpickle\"):\n",
    "        print_color((('To run another Spread, click Reset', 'red'),))\n",
    "        os.remove(\"G_inf.gpickle\")\n",
    "        os.remove(\"data_saver.mat\")\n",
    "        clear_output(wait=True)\n",
    "    else:\n",
    "        print_color((('Click the menu Run Button', 'green'),))\n",
    "\n",
    "buttonres.on_click(reset)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Running the spreading program yields two graphs: a grid of infected squares (left) and relative plots of the perimeter function $\\frac{perim}{4N}$ (blue curve) and infected squares (red curve).\n",
    "[Goto the full spreading code](full.ipynb)"
   ]
  },
  {
   "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.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}