{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# cadCAD Tutorials: The Robot and the Marbles, part 1\n",
    "## What is cadCAD?\n",
    "cadCAD is a Python library that assists in the processes of designing, testing and validating complex systems through simulation. At its core, cadCAD is a differential games engine that supports parameter sweeping and Monte Carlo analyses and can be easily integrated with other scientific computing Python modules and data science workflows.\n",
    "\n",
    "At [BlockScience](http://bit.ly/github_articles_M_1), we've been using cadCAD as a tool for [Token Engineering](https://medium.com/block-science/on-the-practice-of-token-engineering-part-i-c2cc2434e727) - the design of self-organizing systems enabled through cryptographic peer-to-peer networks. But cadCAD can simulate any system that can be described as state variables that evolve over time according to a set of equations. This series of articles will go over the basic concepts of cadCAD and the 'language' in which a system must be described in order for it to be interpretable by the library. This article will cover the following concepts:\n",
    "* State Variables\n",
    "* Timestep\n",
    "* State Update Functions\n",
    "* Partial State Update Blocks\n",
    "* Simulation Configuration Parameters\n",
    "\n",
    "To help illustrate, we'll refer to a simple example system\n",
    "\n",
    "__The robot and the marbles__  \n",
    "* Picture a box (`box_A`) with ten marbles in it; an empty box (`box_B`) next to the first one; and a robot arm capable of taking a marble from any one of the boxes and dropping it into the other one. \n",
    "* The robot is programmed to take one marble at a time from the box containing the largest number of marbles and drop it in the other box. It repeats that process until the boxes contain an equal number of marbles.  \n",
    "\n",
    "## State Variables\n",
    ">A state variable is one of the set of variables that are used to describe the mathematical \"state\" of a dynamical system. Intuitively, the state of a system describes enough about the system to determine its future behaviour in the absence of any external forces affecting the system. [_(source: Wikipedia)_](https://en.wikipedia.org/wiki/State_variable)\n",
    "\n",
    "cadCAD can handle state variables of any Python data type, including custom classes. It is up to the user of cadCAD to determine the state variables needed to sufficiently accurately describe the system they are interested in.\n",
    "\n",
    "We could describe the simple system in our example with only two state variables: the number of marbles in `box_A` and in `box_B`. These are not the only two variables that describe the system, of course. Things like the position of the robot arm in space or its temperature also qualify as \"variables that describe the state of the system\". But if we assume those variables have no impact on the behavior of the system (as implied by the description) we can safely disregard them.\n",
    "\n",
    "cadCAD expects state variables to be passed to it as a Python `dict` where the `dict_keys` are the __names of the \n",
    "variables__ and the `dict_values` are their __initial values__."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
    "# List of all the state variables in the system and their initial values\n",
    "initial_conditions = {\n",
    "    'box_A': 10, # as per the description of the example, box_A starts out with 10 marbles in it\n",
    "    'box_B': 0 # as per the description of the example, box_B starts out empty\n",
    "}\n",
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Timestep\n",
    "Computer simulations run in discrete time:\n",
    ">Discrete time views values of variables as occurring at distinct, separate \"points in time\", or equivalently as being unchanged throughout each non-zero region of time (\"time period\")—that is, time is viewed as a discrete variable. Thus a non-time variable jumps from one value to another as time moves from one time period to the next. This view of time corresponds to a digital clock that gives a fixed reading of 10:37 for a while, and then jumps to a new fixed reading of 10:38, etc. [_(source: Wikipedia)_](https://en.wikipedia.org/wiki/Discrete_time_and_continuous_time)\n",
    "\n",
    "The concept of Timestep in cadCAD refers to a discrete unit of time. cadCAD increments a \"time counter\", and at each step it updates the state variables according to the equations that describe the system. If time itself is a state variable that the user is interested in keeping track of, they may specify a state update function that models the passing of time. We'll cover that scenario in an upcoming article. For the moment, it suffices to define that the robot acts once per timestep.\n",
    "# State Update Functions\n",
    "State Update Functions are Python functions that represent the equations according to which the state variables change over time. Each state update function must return a tuple containing the name of the state variable being updated and its new value. The general structure of a state update function is:\n",
    "```python\n",
    "def function(params, step, sL, s, _input):\n",
    "    ...\n",
    "    y = ...\n",
    "    x = ...\n",
    "    return (y, x)\n",
    "```\n",
    "State update functions can read the current state of the system from argument `s`. We'll ignore the other arguments for now. `s` is a Python `dict` where the `dict_keys` are the __names of the variables__ and the `dict_values` are their __current values__. With this, we can define the state update functions for variables `box_A` and `box_B`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update_A(params, step, sL, s, _input):\n",
    "    y = 'box_A'\n",
    "    add_to_A = 0\n",
    "    if (s['box_A'] > s['box_B']):\n",
    "        add_to_A = -1\n",
    "    elif (s['box_A'] < s['box_B']):\n",
    "        add_to_A = 1\n",
    "    x = s['box_A'] + add_to_A\n",
    "    return (y, x)\n",
    "\n",
    "def update_B(params, step, sL, s, _input):\n",
    "    y = 'box_B'\n",
    "    add_to_B = 0\n",
    "    if (s['box_B'] > s['box_A']):\n",
    "        add_to_B = -1\n",
    "    elif (s['box_B'] < s['box_A']):\n",
    "        add_to_B = 1\n",
    "    x = s['box_B'] + add_to_B\n",
    "    return (y, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Partial State Update Blocks\n",
    "Within a timestep, state update functions can be run in any combination of serial or parallel executions. Take the following diagram for example:\n",
    "![partial state update blocks](partial-state-update-blocks.png \"Partial State Update Blocks\") <center><b>Figure 1: Visual representation of Partial State Update Blocks</b></center>\n",
    "\n",
    "State update functions (SUF) 1 and 2 are run in parallel. This means that if SUF2 reads the value of variable A, it will not get the value updated by SUF1. On the other hand, SUF3 and SUF4 are executed after SUF1 and SUF2 have completed, thus having access to the updated values of variables A and C.\n",
    "\n",
    "We refer to the groups of state update functions that are executed in parallel within a timestep as Partial State Update Blocks. cadCAD expects partial state update blocks to be specified as a list of `dict`s with the following structure:\n",
    "```python\n",
    "partial_state_update_blocks = [\n",
    "    { \n",
    "        'policies': {\n",
    "            'policy1': policy_function_1,\n",
    "            'policy2': policy_function_2,\n",
    "            ...\n",
    "        },\n",
    "        'variables': {\n",
    "            'variable1': state_update_function_1,\n",
    "            'variable2': state_update_function_2,\n",
    "            ...\n",
    "        }\n",
    "    },\n",
    "    ...\n",
    "]\n",
    "```\n",
    "\n",
    "We'll ignore the `policies` key for now. The `dict` that represents the structure of Figure 1 would be:\n",
    "```python\n",
    "partial_state_update_blocks = [\n",
    "    { \n",
    "        'policies': {\n",
    "        },\n",
    "        'variables': {\n",
    "            'variableA': state_update_function_1,\n",
    "            'variableC': state_update_function_2,\n",
    "        }\n",
    "    },\n",
    "    { \n",
    "        'policies': {\n",
    "        },\n",
    "        'variables': {\n",
    "            'variableA': state_update_function_3,\n",
    "            'variableB': state_update_function_4,\n",
    "        }\n",
    "    }\n",
    "]\n",
    "```\n",
    "\n",
    "In the case of our robot and marbles example system, we can model the system so that all state update functions are executed in parallel. In other words, we consider the marbles move from one box to the other simultaneously (ie, `box_A + box_B` is constant)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
    "# In the Partial State Update Blocks, the user specifies if state update functions will be run in series or in parallel\n",
    "partial_state_update_blocks = [\n",
    "    { \n",
    "        'policies': { # We'll ignore policies for now\n",
    "        },\n",
    "        'variables': { # The following state variables will be updated simultaneously\n",
    "            'box_A': update_A,\n",
    "            'box_B': update_B\n",
    "        }\n",
    "    }\n",
    "]\n",
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulation Configuration Parameters\n",
    "Lastly, we define the number of timesteps and the number of Monte Carlo runs of the simulation. These parameters must be passed in a dictionary, in `dict_keys` `T` and `N`, respectively. In our example, we'll run the simulation for 10 timesteps. And because we are dealing with a deterministic system, it makes no sense to have multiple Monte Carlo runs, so we set `N=1`. We'll ignore the `M` key for now and set it to an empty `dict`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
    "# Settings of general simulation parameters, unrelated to the system itself\n",
    "# `T` is a range with the number of discrete units of time the simulation will run for;\n",
    "# `N` is the number of times the simulation will be run (Monte Carlo runs)\n",
    "# In this example, we'll run the simulation once (N=1) and its duration will be of 10 timesteps\n",
    "# We'll cover the `M` key in a future article. For now, let's leave it empty\n",
    "simulation_parameters = {\n",
    "    'T': range(10),\n",
    "    'N': 1,\n",
    "    'M': {}\n",
    "}\n",
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Putting it all together\n",
    "We have defined the state variables of our system and their initial conditions, as well as the state update functions, which have been grouped in a single state update block. We have also specified the parameters of the simulation (number of timesteps and runs). We are now ready to put all those pieces together in a `Configuration` object."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cadCAD.configuration import Configuration\n",
    "\n",
    "# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
    "# The configurations above are then packaged into a `Configuration` object\n",
    "config = Configuration(initial_state=initial_conditions, #dict containing variable names and initial values\n",
    "                       partial_state_update_blocks=partial_state_update_blocks, #dict containing state update functions\n",
    "                       sim_config=simulation_parameters #dict containing simulation parameters\n",
    "                      )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Running the engine\n",
    "We are now ready to run the engine with the configuration defined above. Instantiate an ExecutionMode, an ExecutionContext and an Executor objects, passing the Configuration object to the latter. Then run the `execute()` method of the Executor object, which returns the results of the experiment in the first element of a tuple."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "from cadCAD.engine import ExecutionMode, ExecutionContext, Executor\n",
    "exec_mode = ExecutionMode()\n",
    "exec_context = ExecutionContext(exec_mode.single_proc)\n",
    "executor = Executor(exec_context, [config]) # Pass the configuration object inside an array\n",
    "raw_result, tensor = executor.execute() # The `execute()` method returns a tuple; its first elements contains the raw results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analyzing the results\n",
    "We can now convert the raw results into a DataFrame for analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th>box_A</th>\n",
       "      <th>box_B</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>run</th>\n",
       "      <th>timestep</th>\n",
       "      <th>substep</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th rowspan=\"11\" valign=\"top\">1</th>\n",
       "      <th>0</th>\n",
       "      <th>0</th>\n",
       "      <td>10</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <th>1</th>\n",
       "      <td>9</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <th>1</th>\n",
       "      <td>8</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <th>1</th>\n",
       "      <td>7</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <th>1</th>\n",
       "      <td>6</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                      box_A  box_B\n",
       "run timestep substep              \n",
       "1   0        0           10      0\n",
       "    1        1            9      1\n",
       "    2        1            8      2\n",
       "    3        1            7      3\n",
       "    4        1            6      4\n",
       "    5        1            5      5\n",
       "    6        1            5      5\n",
       "    7        1            5      5\n",
       "    8        1            5      5\n",
       "    9        1            5      5\n",
       "    10       1            5      5"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import pandas as pd\n",
    "df = pd.DataFrame(raw_result)\n",
    "df.set_index(['run', 'timestep', 'substep'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 15376 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "df.plot('timestep', ['box_A', 'box_B'], grid=True, \n",
    "        colormap = 'RdYlGn',\n",
    "        xticks=list(df['timestep'].drop_duplicates()), \n",
    "        yticks=list(range(1+(df['box_A']+df['box_B']).max())));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because the number of marbles in the system is even, it converges to an equilibrium with 5 marbles in each box. Simulating a scenario with an odd number of marbles is as easy as modifying the `initial_condition` of the system, recreating the configuration object and rerunning the simulation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%capture\n",
    "initial_conditions = {\n",
    "    'box_A': 11,\n",
    "    'box_B': 0\n",
    "}\n",
    "config = Configuration(initial_state=initial_conditions, #dict containing variable names and initial values\n",
    "                       partial_state_update_blocks=partial_state_update_blocks, #dict containing state update functions\n",
    "                       sim_config=simulation_parameters #dict containing simulation parameters\n",
    "                      )\n",
    "executor = Executor(exec_context, [config])\n",
    "raw_result, tensor = executor.execute()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th>box_A</th>\n",
       "      <th>box_B</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>run</th>\n",
       "      <th>timestep</th>\n",
       "      <th>substep</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th rowspan=\"11\" valign=\"top\">1</th>\n",
       "      <th>0</th>\n",
       "      <th>0</th>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <th>1</th>\n",
       "      <td>10</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <th>1</th>\n",
       "      <td>9</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <th>1</th>\n",
       "      <td>8</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <th>1</th>\n",
       "      <td>7</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <th>1</th>\n",
       "      <td>6</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <th>1</th>\n",
       "      <td>6</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <th>1</th>\n",
       "      <td>6</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <th>1</th>\n",
       "      <td>5</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                      box_A  box_B\n",
       "run timestep substep              \n",
       "1   0        0           11      0\n",
       "    1        1           10      1\n",
       "    2        1            9      2\n",
       "    3        1            8      3\n",
       "    4        1            7      4\n",
       "    5        1            6      5\n",
       "    6        1            5      6\n",
       "    7        1            6      5\n",
       "    8        1            5      6\n",
       "    9        1            6      5\n",
       "    10       1            5      6"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.DataFrame(raw_result)\n",
    "df.set_index(['run', 'timestep', 'substep'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 20384 characters]\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "df.plot('timestep', ['box_A', 'box_B'], grid=True, \n",
    "        colormap = 'RdYlGn',\n",
    "        xticks=list(df['timestep'].drop_duplicates()), \n",
    "        yticks=list(range(1+(df['box_A']+df['box_B']).max())));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As was to be expected, the system oscilates between 5 and 6 marbles in each box.\n",
    "\n",
    "In the next article of this series we'll cover another base concept in cadCAD: **policies**.\n",
    "\n",
    "---\n",
    "\n",
    "_About BlockScience_  \n",
    "[BlockScience](http://bit.ly/github_articles_M_1) is a research and engineering firm specialized in complex adaptive systems and applying practical methodologies from engineering design, development and testing to projects in emerging technologies such as blockchain. Follow us on [Medium](http://bit.ly/bsci-medium) or [Twitter](http://bit.ly/bsci-twitter) to stay in touch."
   ]
  }
 ],
 "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
