docs/bigidea_TablArray.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "arbitrary-investigation",
   "metadata": {},
   "source": [
    "# the Big Idea of TablArray: math operators for tables of cells\n",
    "\n",
    "In science and engineering we often write formulas which mix scalar, vector, and matrix quantities.\n",
    "\n",
    "But as good programmers, we want our formulas coded to be just as happy to operate on a table of data rather than a point. This means our parameters need to be arrays of scalars, arrays of vectors, or arrays of matrices.\n",
    "\n",
    "Using tradition array coding, broadcasting (and a tool like einsum) can be instrumental in making this happen. But in some systems, code smells arrise especially because some convention is required to define array format.\n",
    "\n",
    "Even worse, if the a physics / engineering lib is generic, how can all applications agree on a convention, or alternately, should libraries impose a convention on applications?\n",
    "\n",
    "TablArray is like a traditional array, with 3 new ideas that answer the call.\n",
    "\n",
    "## Idea 1. Tabular and cellular shapes shall be distinguishable at TablArray instantiation.\n",
    "\n",
    "Otherwise, you see a 3 dimensional array shape; is it a 1d array of matrices, a 2d array of vectors, or a 3d array of scalars?\n",
    "\n",
    "TablArray's first aim is to resolve this ambiguity from the beginning and then upgrade math operators to match.\n",
    "\n",
    "To illustrate, consider a grid describing some scalar t depending on x, y:\n",
    " \n",
    "| | x=0 | x=.5 | x=1 |\n",
    "|---|---|---|---|\n",
    "|y=-1| -15 | -15.25 | -16 |\n",
    "|y=0| 5 | 4.75 | 4 |\n",
    "|y=1| 25 | 24.75 | 24 |\n",
    "|y=2| 45 | 44.75 | 44 |\n",
    "\n",
    "That looks like a 4x3 array. No ambiguity here. But what about a table of vectors v?\n",
    "\n",
    "| | x=0 | x=.5 | x=1 |\n",
    "|---|---|---|---|\n",
    "|y=-1| [0, 1] | [.7, .7] | [1, 0] |\n",
    "|y=0| [0, 1] | [.7, .7] | [1, 0] |\n",
    "|y=1| [0, 1] | [.7, .7] | [1, 0] |\n",
    "|y=2| [0, 1] | [.7, .7] | [1, 0] |\n",
    "\n",
    "And, to make things more interesting, also note that is degenerate along y, so let's simplify. Then v could be represented as a 3x2 array.\n",
    "\n",
    "| x=0 | x=.5 | x=1 |\n",
    "|---|---|---|\n",
    "| [0, 1] | [.7, .7] | [1, 0] |\n",
    "\n",
    "But now we can't write a simple function like <code>t * v</code> anymore. The shape (3, 2) is not broadcast-compatible with (3, 4), so <code>t * v</code> will return an error (as it should). But you can see the intention of <code>t * v</code>.\n",
    "\n",
    "How do we define \"shape\" to distinguish between tabular and cellular structure? The shape of t could be described like (3, 4 |), then the shape of v would be (3, 1 | 2). '|' is a new symbol to mark where tabular shape ends and cellular shape begins.\n",
    "\n",
    "TablArray's primary benefits will derive from this basic idea.\n",
    "\n",
    "## Idea 2. Broadcasting rules then consider tabular **and** cellular shapes.\n",
    "\n",
    "The first, and most powerful, benefit is an upgrade to broadcasting rules. In other words, tables are broadcast, and at the same time, cells are broadcast.\n",
    "\n",
    "In the above example, the degenerate table shapes (3) and (4, 3) will be broadcast together. And at the same time the scalar cells shapes () and (2) will be broadcast.\n",
    "\n",
    "### Code example\n",
    "\n",
    "Let's code the <code>t * v</code> example from earlier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "exotic-champion",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import tablarray as ta\n",
    "x = ta.TablArray(np.linspace(0, 1, 3), 0)\n",
    "# notice 0 after all the linspace()\n",
    "# this tells TablArray that cdim=0, meaning the cells are 0 dim\n",
    "# TablArray always needs to know cdim\n",
    "y = ta.TablArray(np.linspace(-1, 2, 4).reshape(4, 1), 0)\n",
    "# again cdim=0, so cells are scalars\n",
    "t = 20*y + 5 - 1 * x**2\n",
    "v = ta.TablArray([[1, 0], [.7, .7], [0, 1]], 1)\n",
    "# this time cdim=1, so cells are vectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "adaptive-pontiac",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[|0.  0.5 1. |]t(3,)|c()\n",
      "[[|-1.|]\n",
      " [| 0.|]\n",
      " [| 1.|]\n",
      " [| 2.|]]t(4, 1)|c()\n"
     ]
    }
   ],
   "source": [
    "print(x)\n",
    "print(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "functioning-commerce",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[|-15.   -15.25 -16.  |]\n",
      " [|  5.     4.75   4.  |]\n",
      " [| 25.    24.75  24.  |]\n",
      " [| 45.    44.75  44.  |]]t(4, 3)|c()\n",
      "[|[1.  0. ]|\n",
      " |[0.7 0.7]|\n",
      " |[0.  1. ]|]t(3,)|c(2,)\n"
     ]
    }
   ],
   "source": [
    "print(t)\n",
    "print(v)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "thrown-impression",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[|[-15.     -0.   ]|\n",
      "  |[-10.675 -10.675]|\n",
      "  |[ -0.    -16.   ]|]\n",
      "\n",
      " [|[  5.      0.   ]|\n",
      "  |[  3.325   3.325]|\n",
      "  |[  0.      4.   ]|]\n",
      "\n",
      " [|[ 25.      0.   ]|\n",
      "  |[ 17.325  17.325]|\n",
      "  |[  0.     24.   ]|]\n",
      "\n",
      " [|[ 45.      0.   ]|\n",
      "  |[ 31.325  31.325]|\n",
      "  |[  0.     44.   ]|]]t(4, 3)|c(2,)\n"
     ]
    }
   ],
   "source": [
    "print(t * v)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "stupid-custom",
   "metadata": {},
   "source": [
    "Notice that when we print TablArray parameters, '|' is used to show the separation between tabular and cellular structure. Also there is a full shape descriptor like 't(4, 3)|c(2,)' at the tail of the string.\n",
    "\n",
    "### Alternatives\n",
    "\n",
    "Nothing that TablArray can do is 100% unique.\n",
    "\n",
    "Do you want broadcasting to always work like the <code>t * v</code> example? There is a convention you can use that accomplishes the same thing with arrays of arrays. When you start coding your system, decide what the maximum cell dimension is, probably 1 or 2 (vectors or matrices). Once this decision is made you need to be consistent throughout all application and library code. Whenever you instantiate a parameter, pad the cells shape on the left with 1, up to the max cell dimension.\n",
    "\n",
    "For the working example above, the shape of t would need to be (4, 3, 1), while the shape of v would be (4, 1, 1). Here's how that would look in tradition code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "charitable-costa",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[[-15.     -0.   ]\n",
      "  [-10.675 -10.675]\n",
      "  [ -0.    -16.   ]]\n",
      "\n",
      " [[  5.      0.   ]\n",
      "  [  3.325   3.325]\n",
      "  [  0.      4.   ]]\n",
      "\n",
      " [[ 25.      0.   ]\n",
      "  [ 17.325  17.325]\n",
      "  [  0.     24.   ]]\n",
      "\n",
      " [[ 45.      0.   ]\n",
      "  [ 31.325  31.325]\n",
      "  [  0.     44.   ]]]\n"
     ]
    }
   ],
   "source": [
    "x_olde = np.linspace(0, 1, 3).reshape(3, 1)\n",
    "y_olde = np.linspace(-1, 2, 4).reshape(4, 1, 1)\n",
    "t_olde = 20*y_olde + 5 - 1 * x_olde**2\n",
    "v_olde = np.array([[1, 0], [.7, .7], [0, 1]])\n",
    "print(t_olde * v_olde)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "quiet-palace",
   "metadata": {},
   "source": [
    "That's it for broadcasting. Now for the final basic principle of TablArray's:\n",
    "\n",
    "## Idea 3. Matrix operations like matmul always align to cells.\n",
    "\n",
    "In traditional stacks of arrays, you would need to use einsum to implement matmul safely. In fact, tablarray uses numpy.einsum under the hood for matmul. You may have already known about this, and it may not seem like a big deal to do it yourself, even though some engineers will find einsum pretty hard to read. But more often, engineers don't know this and spend time writing outer loops to descend into their stacks and execute matmul on the inner arrays. That's executionally inefficient, as well as inflexible toward varying tabular shape. Experienced coders know to avoid loops and write vectorized code instead.\n",
    "\n",
    "TablArray's table/cell distinction is needed to make matmul always unambiguous for any pair of compatible parameters. Ultimately what this delivers to all engineers is the power of a name. If what you mean is matmul (or the binary operator '@') - call it what it is.\n",
    "\n",
    "## Conclusions\n",
    "\n",
    "I would suggest users first consider the shape of their parameters, and variability in those shapes, to see whether they have a use case for TablArray.\n",
    "\n",
    "TablArray is not functionally unique but may be more readable, bug resistant, and extensible by way of being unambiguous and natural. In other words, traditional stacks of arrays may be abused as often as well used.\n",
    "\n",
    "I would urge the reader to also consider the next levels. TablArray enables more big ideas. Head back to the index to learn about:\n",
    "\n",
    "1. TablaSet, a key-indexed dataset of mutually broadcast-compatible parameters. It bears resemblance to a database, but with fast advanced slicing called projection.\n",
    "2. TablaSolve, an object oriented dataset where user-defined solver methods are attached at run-time. Anytime seed data is presented or manipulated, solvers can be automatically arranged to determine downstream parameters.\n",
    "\n",
    "## More info: indexing\n",
    "\n",
    "TablArray objects have a property .view which determines the behavior of indexing. Views make use of pass-by-reference. For example if t1 is a TablArray then\n",
    "\n",
    "* t1.cell[0] will return a view of t1's full table looking at all the 0th elements of the cells.\n",
    "* t1.table[0] will return a view of t1's full cell at the 0th table position.\n",
    "* t1.bcast is like table view but uses projection - there is a better discussion under TablaSet. Technically what this means is that if you index missing dimensions or dimensions which have 1 length, you will get an answer instead of an error. The answer you get will depend on broadcasting rules during indexing. But the term 'projection' makes more tactile sense: if you call for a row of a 2d table, you get a row, while if you call for a row of a 1d table, you get a scalar. That's not the same as indexing; where if you call for a row of a 1d table you get an error.\n",
    "* t1.array is a view of t1 which will ignore tabular and cellular distinction. E.g. a shape like t(4, 3)|c(1) will be indexed the same as a traditionally (4, 3, 1) array.\n",
    "\n",
    "Example code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "simple-preserve",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.7 0.7]\n"
     ]
    }
   ],
   "source": [
    "print(v.table[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "compound-actress",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[|0.  0.7 1. |]t(3,)|c()\n"
     ]
    }
   ],
   "source": [
    "print(v.cell[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "wooden-macro",
   "metadata": {},
   "source": [
    "TablArray objects also have a .base property which is actually the numpy.ndarray data that implements the TablArray object. In case you find the need to mix TablArray and other ndarray methods."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "env1",
   "language": "python",
   "name": "env1"
  },
  "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}