2017-03-17_mg_bs_driven_feedback_setup.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from collections import OrderedDict\n",
    "from copy import copy\n",
    "from sympy import Symbol\n",
    "import QDYN\n",
    "from QDYN.pulse import blackman\n",
    "from QDYN.shutil import mkdir\n",
    "import matplotlib.pylab as plt\n",
    "import numpy as np\n",
    "import sympy\n",
    "from sympy import pi\n",
    "\n",
    "from src import single_sided_network_v1 as single_sided_network\n",
    "from src import qdyn_model_v1 as qdyn_model\n",
    "from src.notebook_plots_v1 import display_hamiltonian"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Circuit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "circuit = single_sided_network.network_circuit(n_nodes=2, topology='driven_bs_fb')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAscAAABhCAIAAAC1TpujAAAJNmlDQ1BkZWZhdWx0X3JnYi5pY2MA\nAHiclZFnUJSHFobP933bCwvssnRYepMqZQHpvUmvogJL7yxLEbEhYgQiiog0RZCggAGjUiRWRLEQ\nFBSxoFkkCCgxGEVUUPLDOxPn3vHHfX49884755yZA0ARBQBARQFSUgV8Pxd7TkhoGAe+IZKXmW7n\n4+MJ3+X9KCAAAPdWfb/zXSjRMZk8AFgGgHxeOl8AgOQCgGaOIF0AgBwFAFZUUroAADkLACx+SGgY\nAHIDAFhxX30cAFhRX30eAFj8AD8HABQHQKLFfeNR3/h/9gIAKNvxBQmxMbkc/7RYQU4kP4aT6edi\nz3FzcOD48NNiE5Jjvjn4/yp/B0FMrgAAwCEtfRM/IS5ewPmfoUYGhobw7y/e+gICAAh78L//AwDf\n9NIaAbgLANi+f7OoaoDuXQBSj//NVI8CMAoBuu7wsvjZXzMcAAAeKMAAFkiDAqiAJuiCEZiBJdiC\nE7iDNwRAKGwAHsRDCvAhB/JhBxRBCeyDg1AD9dAELdAOp6EbzsMVuA634S6MwhMQwhS8gnl4D0sI\nghAROsJEpBFFRA3RQYwQLmKNOCGeiB8SikQgcUgqkoXkIzuREqQcqUEakBbkF+QccgW5iQwjj5AJ\nZBb5G/mEYigNZaHyqDqqj3JRO9QDDUDXo3FoBpqHFqJ70Sq0ET2JdqFX0NvoKCpEX6ELGGBUjI0p\nYboYF3PAvLEwLBbjY1uxYqwSa8TasV5sALuHCbE57COOgGPiODhdnCXOFReI4+EycFtxpbga3Alc\nF64fdw83gZvHfcHT8XJ4HbwF3g0fgo/D5+CL8JX4Znwn/hp+FD+Ff08gENgEDYIZwZUQSkgkbCaU\nEg4TOgiXCcOEScICkUiUJuoQrYjexEiigFhErCaeJF4ijhCniB9IVJIiyYjkTAojpZIKSJWkVtJF\n0ghpmrREFiWrkS3I3uRo8iZyGbmJ3Eu+Q54iL1HEKBoUK0oAJZGyg1JFaadco4xT3lKpVGWqOdWX\nmkDdTq2inqLeoE5QP9LEado0B1o4LYu2l3acdpn2iPaWTqer023pYXQBfS+9hX6V/oz+QYQpoifi\nJhItsk2kVqRLZETkNYPMUGPYMTYw8hiVjDOMO4w5UbKouqiDaKToVtFa0XOiY6ILYkwxQzFvsRSx\nUrFWsZtiM+JEcXVxJ/Fo8ULxY+JXxSeZGFOF6cDkMXcym5jXmFMsAkuD5cZKZJWwfmYNseYlxCWM\nJYIkciVqJS5ICNkYW53txk5ml7FPsx+wP0nKS9pJxkjukWyXHJFclJKVspWKkSqW6pAalfokzZF2\nkk6S3i/dLf1UBiejLeMrkyNzROaazJwsS9ZSlidbLHta9rEcKqct5ye3We6Y3KDcgryCvIt8uny1\n/FX5OQW2gq1CokKFwkWFWUWmorVigmKF4iXFlxwJjh0nmVPF6efMK8kpuSplKTUoDSktKWsoByoX\nKHcoP1WhqHBVYlUqVPpU5lUVVb1U81XbVB+rkdW4avFqh9QG1BbVNdSD1Xerd6vPaEhpuGnkabRp\njGvSNW00MzQbNe9rEbS4Wklah7XuaqPaJtrx2rXad3RQHVOdBJ3DOsOr8KvMV6Wualw1pkvTtdPN\n1m3TndBj63nqFeh1673WV9UP09+vP6D/xcDEINmgyeCJobihu2GBYa/h30baRjyjWqP7q+mrnVdv\nW92z+o2xjnGM8RHjhyZMEy+T3SZ9Jp9NzUz5pu2ms2aqZhFmdWZjXBbXh1vKvWGON7c332Z+3vyj\nhamFwOK0xV+WupZJlq2WM2s01sSsaVozaaVsFWnVYCW05lhHWB+1Ftoo2UTaNNo8t1WxjbZttp22\n07JLtDtp99rewJ5v32m/6GDhsMXhsiPm6OJY7DjkJO4U6FTj9MxZ2TnOuc153sXEZbPLZVe8q4fr\nftcxN3k3nluL27y7mfsW934Pmoe/R43Hc09tT75nrxfq5e51wGt8rdra1LXd3uDt5n3A+6mPhk+G\nz6++BF8f31rfF36Gfvl+A/5M/43+rf7vA+wDygKeBGoGZgX2BTGCwoNaghaDHYPLg4Uh+iFbQm6H\nyoQmhPaEEcOCwprDFtY5rTu4bircJLwo/MF6jfW5629ukNmQvOHCRsbGyI1nIvARwRGtEcuR3pGN\nkQtRblF1UfM8B94h3qto2+iK6NkYq5jymOlYq9jy2Jk4q7gDcbPxNvGV8XMJDgk1CW8SXRPrExeT\nvJOOJ60kByd3pJBSIlLOpYqnJqX2pymk5aYNp+ukF6ULMywyDmbM8z34zZlI5vrMHgFLkC4YzNLM\n2pU1kW2dXZv9ISco50yuWG5q7uAm7U17Nk3nOef9tBm3mbe5L18pf0f+xBa7LQ1bka1RW/u2qWwr\n3Da13WX7iR2UHUk7fiswKCgveLczeGdvoXzh9sLJXS672opEivhFY7std9f/gPsh4YehPav3VO/5\nUhxdfKvEoKSyZLmUV3rrR8Mfq35c2Ru7d6jMtOzIPsK+1H0P9tvsP1EuVp5XPnnA60BXBaeiuOLd\nwY0Hb1YaV9YfohzKOiSs8qzqqVat3le9XBNfM1prX9tRJ1e3p27xcPThkSO2R9rr5etL6j8dTTj6\nsMGloatRvbHyGOFY9rEXTUFNAz9xf2pplmkuaf58PPW48ITfif4Ws5aWVrnWsja0Latt9mT4ybs/\nO/7c067b3tDB7ig5BaeyTr38JeKXB6c9Tved4Z5pP6t2tq6T2VnchXRt6prvju8W9oT2DJ9zP9fX\na9nb+aver8fPK52vvSBxoewi5WLhxZVLeZcWLqdfnrsSd2Wyb2Pfk6shV+/3+/YPXfO4duO68/Wr\nA3YDl25Y3Th/0+LmuVvcW923TW93DZoMdv5m8lvnkOlQ1x2zOz13ze/2Dq8ZvjhiM3LlnuO96/fd\n7t8eXTs6/CDwwcOx8DHhw+iHM4+SH715nP146cn2cfx48VPRp5XP5J41/q71e4fQVHhhwnFi8Ln/\n8yeTvMlXf2T+sTxV+IL+onJacbplxmjm/Kzz7N2X615OvUp/tTRX9KfYn3WvNV+f/cv2r8H5kPmp\nN/w3K3+XvpV+e/yd8bu+BZ+FZ+9T3i8tFn+Q/nDiI/fjwKfgT9NLOcvE5arPWp97v3h8GV9JWVn5\nBy6ikLxSF1/9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAAHXRFWHRTb2Z0d2FyZQBHUEwgR2hvc3Rz\nY3JpcHQgOS4xOJQFEHMAABc0SURBVHic7d1/jBxneQfwx+FwYuPYfq+gOqSKj1kTF5EfdubsRi3U\n12imAU6tAmZXDkQgRJhpVKmqVKk7ImplCRVm1FKgVUV3UgtVGCXsKEpEeqTmXhkbSmTwTWKnBkGc\nHXxRqgtpfXP+IYMtJ+4frzPeu513bvc8uzvv7vfzR/Te3O7m8c3u7DPv+7zvu+rKlSskE4YUx6Tr\nxBhxTkTStmEQ0bV2HFMYLmozRrq+qB1FFEXStqaRpi1qi2Ca2yKAwQuyNbDrD7I1sE6DXPav12mQ\n7fz12g7yl1/Yf/nVuct/eN/7Pqm//o/7z/xi7sZP7L7tLvb6N6aS9uzfP3Hp9Pmk/esLtOkv92zc\nQNfab8YnvSffWLMuaV8cveV37MmkfceXHqIwPOkfFO1z/3HotWeONrffHN+x9XMTzW0RTHN7JUG+\nGc/+24H2g2wN7MKt7737bx9oDuz8t56e+8FJWZA3fmL3bbu05iBbA1v70fs7CvKOLz1EUXT6sSd/\n6/230EMP0aFDdPQo3dLUfu976YEHrrZ37KCJCXr6aTp58mp7/36am1vU3r2bNO1amzF64gk6f/5a\nm4juv/9ae88eimM6cOBa+8knad26a20RTBRdaxc2yDCkgwc7CLI1sPvuI11fFNjUVFaQe/YQ0aIg\nWwObnGwzyFefOvrabTtGPzqhvXLo9IGjv3znjtGPTmjP7X/16NzV48/tf/Xo3P/87n3vul8X7dkd\nu9//Abbxv6autZ994uSx869/8Gr757+g1yb2TOwieuKt9l3xhf1PvjK/LmlHv77l/35vMmlf/PhD\nOoUXpg5ebZ87dPrA0Tm61r4aWOdBnvzhXBJYSpB3xXTgQEaQl++f3HZbU5AtgaUE+eLTp398MiPI\n2/fobGp/xl/v9Pj9rUG+NnLrps/9qfQCnuYGaUoBoI5Xn3m+3yFAW378435HAAXw4n/3OwJoz6v/\n/FT6L8RNXZpVWX0VnF+9QwUotoVHnI1fd/sdBSzP98my+h0E9Ju414Xie/yL0YOf7+xUZfZVyJMR\ngELZuBvprxqQUgARUgplSFMKeXowkvV6aWc+juNQDKJDz2mapuHjmMbhhpuWV3Bkxn3CGNN1vfV4\nLjepURRFUXS9rwIrous6Y+w6X8RxyEXfogqk3cDyoYzMrCLttiIMQ8/zUq8X0FVxHDPGXHwW02yZ\n5UQpb3Hbtsvlcu/jgTAMp6enW4/7fg5fJ77vi4/D9b4QdCgMw2q1amBkfGgcO0YTqb+QvwcyswpJ\nMqLrOr7beo9zjjtvmYc3p2cVmqbhvdoXpmmmHs+rr61cLuO7rfccx8nldfChVMVLn3EnUn8h//Sh\nrgIGAr5gFIG6CgCFSD+w8vQgM6vAlRoU4XC8V9WAWxUgopy6PKDrTj0oOVXIKmCw7Tzu9zsEaAuy\nCgCFnDol+QVGQJYIgmDZGnLP8+I47k08cJ0+djdmBKgBc5iAUFehjpc+IzlVA5lVcM5N0xwdHRUF\nRFEU2ba9atUqx3FExsA5Hx0dtW17yRPDMIzjeNkpmpZl+T7ugNUwdwf61dSAugogIswLVoVlSE6V\n/MtR4REQwzAsy0ompGiaJuYQuq4rMgbDMFzXrdVqS57o+77VxrVNzFvDtHglfO1Eod+rjuOYpmma\nZtL75fv+6OhopVJpZ16PeHCXY+yRYt+qQI/gfk0Vpx6VnCr5N2NWVhEyVvAhAMZY65JczRfu1uyB\nc97+YhvlchndFUq488jXipz/ua6r63q1Wk0m5lmWZVlWvV5vZ3qkyJ67HGOPdJpVYDY1QB/lXFex\n9lvfKpVKyYBCMSU5RBiG4sor8owoilIvxK1ZhVjXi3PuOI54bvKamqZhIVElfPLeOdM0K5VKEAT9\njkXKMAzGWJEj7IFO6yo458W/CkGnUFehipzrKt63erWu657nlUol0zQLeNfenB+ItfaSaglZn0SS\nfCQ/ep4nVouzLEsUYeD2SDlzdxjlcjkIgkqlUiqVCltpW61WPc9r/YL0PC8IgiAImiMPw9BxHM65\nWEoyebDjOL7vB0HgOE4x/5kZOq2rqFarcRyLq5Bt2/hsDgakiKpYQV3FCBH5vp96H7BldnbdunWi\nLRZ29Dxv586dmzZtyiHYPCSlD1EUiZ5kTdPETiVtrtPseV7SBS16JoIgKOyafYcPH25/bbszZ85s\n2LBBLNpf2H9Rp77zne8899xzrccPH/7QvfdevPnmm8+dOxdFkeM4juNYlnXu3LneB5lBrLlu23bz\nataO4yRVQVEUmaY5MzNDRLZt12o1kQQnRcfiDSCG9jRNq1QqqQtj9938/Hzqe3V2dsvmzS+3Hl+/\nfv3Zs2dTX2rr1q1HjhwhIt/3fd/XNO32228fmLe0Ws6cObNv3772c7tLly6tXr269fjhwx/ates/\ncw0NuuL3n3nxuT+5q/X4H33ve99PSxtc1x0hIl3XZRMi3vntbzf/aBjGli1bTp8+nUe0ueGcN19i\nxJhFm+PQQRC0lnMWdn+BsbGx9i+m58+fP3LkCOdcfFGJncnEudY0TdGL8le+8pVHH3209fg94Q8u\n6fpXv/rV5Iiu67quF7Db3DAMMdyWHPF9v16vi7ZIbcU0pSX9asmDq9VqclkXDfFf8Y9tpxK5B9as\nWZP6Htu3b4thjLUeHxkZuXz5cupLnThxImmLFPmGGzLLzKFrbrrppnvuuWf79u1tPv7ChQtr165t\nPX7p0m8regkaNr88uMEwdrQef/M97zFKpdSnXM0q0l/P9+2DB+mt6RXVapUxVrTdKHRdb54myhiL\noihj6wdd15sv1pqmJTmEuCgXeTOqzZs3d/RRfOCBB5K2+JYKwzCKIjFSwBirVqsF+RJq08jISPpf\ngHP7hz8UTVEIKU5xMYsYqtXq+Ph48iZsHcXIHtfQdT35I1y5coWIHMcR3RumaRbkhMqyCiJKzSpk\noiiampoiIk3TqtVquVxmjOW1GwV06sYbb9y+ffv1JwSGQUTvyyMi6C5OkrO9srqKY089pWlavV5v\nNBqu6xbzDl7X9ebLqLj0ZDxe3CkmP9ZqNVGqKYaBLMsKgiCp0JSVfKpI3ORVq1XXdaenp+fn52u1\nWhAE4+PjhUoTV+ZnH3wwDEPXdcW/q+BnjTEm/vjix+a+hyAIRNJgGIamaa19LUuWUREvIlIKJWos\nOv0+8jzPsqzp6elGo2FZVjGvQtCp4vUhQjpDW1FdhUzJtqebbnaLacn4xbIbVBqG0Xz/Ki7fzT82\nPzgIgoLc+XWD+Ldzzm3bFtmYun2Sj/Gt4pu1mER9pehFSxZTSdJf13VFIQgRxXGc1EnU63Xf9w3D\nEOmC4zhiCRbHcWzbZowtqZhpHkkpLN/vrGCzdYASBoDvYxqIGo79ub/t2bRTJU8Ms7KKmwufUqyM\nuNtbNl0Ql/Jll+BUnWEYjUbD933btsVXnYr/5E//fC9Rca9Sruu25rvNR1KzYVEaItpipCPjwZxz\ny7I450UewiPcpAIRYeF2dUSkbUv9hfwLNLPoqXhTSXMhpkUsW8fX5hKcg8GyrEajwRgbHx9Xcb7i\ntrv7HUFfmabpeV6bK3X2l7LdYZCnobmyKm/9X0lOlTwxzMwqBve2olwuL3tHLqpTexNPQbiu22g0\n4jhuXltaCT+6Y6ivUtNvKf54AbIKICzcrg6DJKdKXjGdmVUgnxw+opBQ13W1EotnTqBHVQ0D2gEK\nnUFWoYpjX+74VGVmFRj7Gla1Ws0wDIUSiz3HMdVQDYPbAQodwHeLKiKSnKoV1lVgUvgQc13XsqxS\nqaTETihDXlehEIyAAKEfXB1511XAcLMsy3Vd0zSLn1gMeV2FQpBVAGEERB0rqKvImlkqyyeb14mC\nnonjuPfrSYhZMKZpTk9PF3lpqWdOaH+QdjwMQ9M0ex0NyHW6XgUMJM6RX6rh2Jf5tg5PVWZWkdbF\nIZY36Oj/AUoTq19XKpV6vV7YxGLPcSd1vYr5+fneBwMZUFcBhLoKdeS9XgXqKoCIiHRdr1arlUql\nsMWbqKtQBToqgPA2UMe2f0FdBXSHZVnlcrmwowncKO7CmtAMN6lAqKtQhxb1pK4ChpPrurZt27Zd\nwKWWMEyriqGqq+CcB0EgNq8XR8SGyc1rrjeXqYlGshHMAMMHVhW9qKuAYSamhHiel70xbO/92Wx6\nXQUUzVDVVYid4cSucuJIHMejo6PlclmUKDmOwxhLfhsEwbI7CQwGpBSqmNlopNdVyHeHy8wqHAf7\nykEzxtj09HSpVNI0rVC7WI1t7ncE0J7h6agQwjBsnroVhiFjLNkuwPf95oJiTdMKWxCdL2QVqrjv\n7zo+VairgM6IxMJxnELNLkZdhSqGrQM0CIIkqwiCwHGcer2ebDDEGGveE07XdRV3DF4B1FWoYgV1\nFZlZBToqIE0Bp4Q8/+RQ9BsPgKGaWCaGM+I45pxzzsMwdF23uetiZmZG0zTP81atWqXQAvnXD1mF\nKg79TcenKnMEBEDCsqwoikzTnJmZ6XcsRESVBR91FVA0nHPDMJI0Qtf10dHRRqMhOiTiOG4uqvB9\n33GcApZCdwNGQFTx8pgxkfoLeacD1quAFXJdV9d127b7HQgR6irUMVR1FZzz5joJ0RWRdEj4i/dv\ntSwrGRkZeMgqVCGtq5CXFaOuAlbOdd0wDD3P63cg5GvoqFDDcJQNXCX6KpIfPc8zDCPJM4IgaB7y\ncBzHGpqca3FCBcWlccmpkp/CzBEQ1FVApuJMCVl4ISLZjr1QJEMyscz3/SiK4jhOOiTEShXJWhTi\nRzGzlIgYY+VyeUhKNWnIJhgr7dA3ookOc93MrCKKhuvOAjonEotKpdLfSXGoq4BCEb0OrjyB0nV9\nSEooUmEERBV511WglwraUIQpIWN/hPRXDcPQUQHLQlahiocfR10F9EOyS0i/EgufhmVAGmAA4I5V\nGbJThboK6DaxS4hpmtPT070vZUddhSqGpK4CssludEdHR4dkddECSt2ARlpXIa+OQF0F5KZWqzmO\n05fEwiLUVQAoQzbZpbmgFXpJth/1T7ZZE6m/kM9XQl0F5Ml1XcuySqVSj9fz3rgd6a8a0FEBNGQT\njJX211+XnKoV1lXgzEPnLMsSW5v2MrFAXYUqMKUQCHes6vjfL3ZcV5GZVQzNkiyQL5FYVCqVniUW\nN3wf+wqoAV8nQEgu1fHTZySnSt7pkJlV4MzDSlmWVa/XK5VKEAQ9+N99fCOyCgBl4I5VFT/ZJjlV\nqKuA3tN1vV6vi/pN3uU9ClFXoQrUVQBhdF0d0roK+SUddRXQRbquNxqNcrls27Zt21HXer+cCPc+\nakAHKFDhd64U90LNC/D4vj86OlqpVNq5QRIP7nKMPbLwiORUrTCrQC8V5MGyrEajwRgbHx+3bbsb\nK2VtmcUIiBrQAQrFJzZkrlarzlvpj2VZYlTXaGNZUMuyBmbVjWPHJL9YYV1Fl3utYai4rityi1Kp\n5DhOvrkF6ioGleM4juN0r5cL+uLDH37FNE2/2DmmYRiMsd5UhhVW3nUVyCogV4wx13VnZmbiOC6V\nSqZpOo7jeR7n/Dq/NlBXoYpO6yosy/I8r1Qq9azyF3pg167bNE2zbVvcYxQ2a6xWq57ntYbneV4Q\nBEEQeJ6X3CCFYeg4Dufc9/0le9z7vh8EQe53Uz2wgrqKzLU1UVcBXaBpWq1Wq9VqIpmIokh8bqMo\nMgxD0zTG2L333rtu3brW5y4sLKSOa+77/h9/VkMSrIDXXluzadOvW4+PjIxcvnw59SmTk5NTU1Pi\nIq5pWrlcvnjxYpfDhHQXL1584YUX2n/8hQsX1q5d23p8374tk5MfFPvFe57neZ7YSyi/SPMhboRs\n225e8dNxHPE+JKIoikzTnJmZISLbtmu1mhj7sG07eTC9tYetpmmVSqWYi4eeO3cu9dJ6yz/9+9xf\nfLr1+Aempm6SDAZlZhWoq4Buah2hDMMwiqIwDB977LE777yz9SmyrOJdL77I+etdiRJydfjwh3bt\nSjmD69evP3v2bOpT3njjjaQdRZHv+1u3bp2cnOxWiCD3m9/85vnnnz99+nSbj7906dLq1atbj586\nNbJkMRuRYczPz+cQZa4Mw+CcO031pb7v1+t10dY0LQzDMAzjOA7DsLWcwvf9arWaXLVEQ7ygyEUK\n4uzZs6mX1vc//+ZPW47/6le/urSw8BHJS2VmFZxjw1roJV3XdV3PuGUJw9BN60NfOONsxJxFFfg+\nWZKNBWRKpZJo6Louiuacgk8hGFwbNmzYvXt3OxWL2aKIOH9JtBlj4rRqmibbjaK/qtXq+Ph4kjG0\njmJkj2voup78xa5cuUJvVWx0IdKVu/XWW1MvrY+vj9zPLx214JxzzmVZBeoqYBBs3I30Vw2ddoCK\nIWrLsmZmZmZmZix0oA4ETSPP83Rdr9VqjUbDdV2twAPujLFarZaU9TT3PQRBIJIGMXrbWoFhWVZz\nUapytUEPtqQUwpbZWdlTUFcBg8Dhhou8QgWdboTMGBNTh7oWEfTBI48s1Gq16+/z6BJRXxlFUZLu\nGIZRrVbFb13XFVOTiCiO46ROol6v+75vGIbounAcxzAM8WDbthljjLHC/pNlFh5xNn49pQ9jy6lT\nsqegrgIGwZZZTqTYx3U4+X5n00AKWMEH12/Dho1F/n51Xbd1OKD5SOpggRjAFW0x0pHxYFUcO0YT\nacdfHhtLPU4YAYHB8PBmvFfVgA5QoKFfuJ1zHsfx9c+o74GXPpN+ql7evFn2FFRrwkDAG1UR6AAF\nMAyjUBNAMsg+sDuPH5c9JbOvAldqUITD8V5VAzpAgQq/DwgkTj2YfqpGFxZkT0FWAYNg5/FCL/0L\nCWQVAAqRFWW+PDYmewrqKmAQfOzuog9PgoC6CqChr6tQyArqKpBVwCCYuwP9ampAXQUQUeGLFOEq\ny0g/VRl1FZnVmhgBAUV87YSBmx8l5FUCvmS9Z+iNvDbHkk0wFjMjcvlfQC5OPeqPPZ5yqjLqKpBV\nwCDYedwnwl2wAnLJKsTWDPj66T3GWFcXwRRntnuvDzKyhWFOnaKxtOMZ61VkZhW+j/5KUALqKlSR\ny62KWCA5hxeCPpHVVSi9YNRAuvyFfOsqMPYFinjlw0h/1YBkAAjfLeowtI7rKt62d+/erJdExTYU\nyTe/+c1PfepTrce/9K9M04gxCkP62c+I6FqbMVqzZlGbc4oiaVvTKI7pRz+StuOY3v3uRe0oojBc\n1BYBNLdbA+s0SPFZHIwgYch5HmV/YJf9kFLm26yd93+nn4UVXEkGIMgb/sHb+PGltwJRFL3j2WfH\nHn449eRivQoYBI+MYL0KNSClACLC7i6qeMed6Z/Yn2zbJnvKquZNUJZCXQUUjGmayQ6BizgOpsAD\nAPSAqJWWFcGgrgIGAtJfAIDcSdKAla5XgSs1FI+DLQQAAHrj8GHatWvJsTiOH5SvV5E5AgJQMGEY\npq/Ds28fffazPQ8HAGCgffe79JGPtB6+/e1vv60l2xBQVwEDAXUVAAAFgLoKGAhIfwEAcidb6lQ+\nEp2ZVeBKDarAhEUAgNx1voB6ZlaBKzWoAiWcAAC5k6UB8k6HzKwCV2oAAIChJcse5J0OmVkFgCow\nWgcAkDvUVcCQwmgdAEDuUFcBQwqjdQAAuZOlAfKZ/KirAAAAgDSdD1mgrgIGApbAAgDIXed1FW/b\nu3ev9OU0jRi73pgAeiCKyPOIczKMlHYUka4T5+T70jYRadqitu9TECxqi09Ec9txiHNp2zCI6Fpb\nFhiC7FeQrYF1GuSygXUaZDt/PQQ5YO/DIgeZvOASsuPYBwQAAADygroKAAAA6IR8Qw/UVQAAAEAn\nxJBNmsysAhVwAAAA0DbUVQAAAEA+UFcBAAAAnUBdBQAAAORDXleRuV6FYRR3prIS06mVCFL16dQI\nEqs4IMiBfB8iyCJ/ouOYdH1pzkBEqKsAAACAvPw/PjhSHIZiFVkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "circuit.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "slh = single_sided_network.network_slh(n_cavity=2, n_nodes=2, topology='driven_bs_fb')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\begin{pmatrix}1\\end{pmatrix}$"
      ],
      "text/plain": [
       "[[1]]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "slh.S"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\begin{pmatrix}- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})}\\end{pmatrix}$"
      ],
      "text/plain": [
       "[[-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁) - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)]]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "slh.L"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def plot_bs_decay():\n",
    "    from qnet.algebra.pattern_matching import pattern, wc\n",
    "    from qnet.algebra.operator_algebra import ScalarTimesOperator\n",
    "    θ = Symbol('theta', real=True)\n",
    "    κ = Symbol('kappa', positive=True)\n",
    "    fig = plt.figure()\n",
    "    ax = fig.add_subplot(111)\n",
    "    theta = np.linspace(0, 2*np.pi, 100)\n",
    "    L = slh.L[0, 0]\n",
    "    k1_form, k2_form = [op.coeff for op in pattern(ScalarTimesOperator).findall(L)]\n",
    "    k1_form = k1_form / sympy.sqrt(2*κ)\n",
    "    k2_form = k2_form / sympy.sqrt(2*κ)\n",
    "    k1 = np.array([k1_form.subs({θ: val}) for val in theta])\n",
    "    k2 = np.array([k2_form.subs({θ: val}) for val in theta])\n",
    "    ax.plot(theta/np.pi, k1)\n",
    "    ax.plot(theta/np.pi, k2)\n",
    "    ax.set_ylim(-3, 3)\n",
    "    ax.set_xlabel(r'BS mixing angle $\\theta$ ($\\pi$ rad)')\n",
    "    ax.set_ylabel(r'decay rate relative to open network')\n",
    "    plt.show(fig)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAh4AAAF9CAYAAABLS2EjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3Xd4VNXaxuHfm0AIBEiQIjbEfhCPNBufBRQrKHYkFhQV\nwS4WEEERLAiCSLGgSBMIIoqgVAVFxYKCigqCohQLIFWpIZn1/THDMWCAZGfP7Mzkua8rFzNr78l+\nPHOSebP2KuacQ0RERCQWkoIOICIiIiWHCg8RERGJGRUeIiIiEjMqPERERCRmVHiIiIhIzKjwEBER\nkZhR4SEiIiIxo8JDREREYkaFh4iIiMSMCg8RERGJmbgrPMysnZl9Y2YbI1+fmNn5QecSERGRfbN4\n26vFzJoBucBPkaYbgAeAus65hUHlEhERkX2Lu8IjP2a2FrjfOTc06CwiIiKyZ6WCDlAUZpYEtADK\nAZ8GHEdERET2IS4LDzM7jnChkQr8DVzqnPsh2FQiIiKyL3F5q8XMSgE1gAzgcqANcEZ+xYeZVQbO\nA5YC22IYU0REJN6lAjWBac65tX58w7gsPHZnZu8CPznnbs3n2NXAqNinEhERSRjXOOdG+/GN4vJW\nSz6SgDJ7OLYUYOTIkdSqVStmgYqrU/tcRzoHMfm+p4KOUign9m3G8WWb8Uq722jfvj19+/YNOpL4\nSO9pYtH7mTgWLlzItddeC5HPUj/EXeFhZk8AU4AVQAXgGqARcO4eXrINoFatWtSvXz8mGYur735Z\nxbajF3BdlY5x97+F7V+G8mkHUL9+fdLT0+Muv+yd3tPEovczIfk2VCHuCg9gf2AEcACwEZgPnOuc\nmxloqjjQb/JUcEb7iy4IOoqIiJRQcVd4OOduDjpDvJq6ZBJp7iRq1agadBQRESmh4m7JdPFmy7Yd\n/Jo6jZMymgYdRURESjAVHiXEy9M+gTJ/ceNpzYKOUmSZmZlBRxCf6T1NLHo/ZW/i7laLeDNqziSS\ndlSnZeN6QUfxbOfUb/1SSzx6TxOL3k/ZG/V4lBDzt07iCNeUUsnx+ZYbFnQEERHxQXx+CkmhTPjk\ne7anL+DK4y8OOoqIiJRwutVSAjw9NQvbkUHHy84LOoqIiJRwKjwSXCjk+HxLFkfbFVRM29PiriIi\nIrGhWy0Jbui7c8ip8DNtGmqwl4iIBE89Hglu4AdZJOUewJ0XNQo6ioiIiHo8Eln2jlzm575GneSr\nSCmdHHScInPE/07KIiIlnXo8Eli/iR8QSlvJnQlxm0XTaUVEEoEKjwT28mejKRU6guvPPjHoKCIi\nIoButSSsvzZv58dSb3BK+UySktRbICIixYN6PBJUj3FTIHUjD5yXCLdZREQkUajwSFCjvski1dWh\n+SnHBh1FRETkf3SrJQH9vvZvVpSbSOMq6u0QEZHiRYVHAnps7AQovY0ul7QMOoqvdu5OKyIi8Uu3\nWhLQ2EXDqcCpnFr70KCj+EgDZEVEEoEKjwTz7twfWVfpPdpUGR50FBERkX9R4ZFgOr3xIhbaj17X\ntwg6ioiIyL9ojEcCWffXVua5odRPbk1G+dSg44iIiPyLejwSSMdXx+JS19Ozebugo4iIiORLhUcC\nyfrpBSq7c2lS78igo4iIiORLhUeCGDVzHpszPufOmuODjhI12p1WRCT+qfBIEI9NeYFkDqZr5oVB\nR4kOp+m0IiKJQINLE8Dy1RtZlDKaM8q3ITVFtaSIiBRfhS48zGyP63Cb2dNFiyNe3Dt8BCRn0/vq\nm4OOIiIisldeejxeNLOmuzeaWV/g2qJHksIIhRxv//ECB226hPpHHRh0HBERkb3yUnhcDYw2s9N3\nNpjZAKAFcKZfwaRg+k+cRXb6QtqfdmvQUURERPap0AMCnHOTzOw2YIKZnQPcBFwMnOmcW+x3QNm7\nx2c9Rar9l/aXqOYTEZHiz9NIROfcaDPLAGYDfwKNnHM/+ZpM9mnkjLmszZjGnQdkkZSkWR8iIlL8\nFajwMLNn9nBoNfAVcJtZ+IPPOXevP9FkXx6c1IPSdiS9H7oy6Cgx4ZzW8RARiXcF7fGot4f2JUDF\nPMf1yRAj73y+kN8qvkmrSi+RUjo56DhRZ6hHR0QkERSo8HDOaQBBMXPP6z1JsgMZ8MB1QUcREREp\nsELNajGzUmaWY2bHRSuQ7NvH3y1lSdpImle5n4ppZYKOIyIiUmCFKjyccznAciDx+/aLsdtG9sa2\nZzDoljZBRxERESkUL+t4PAE8aWb7+R1G9m3+zyv5tvRgzip3D9UqpQUdR0REpFC8TKe9AzgS+N3M\nlgGb8x50ztX3I5jkr82QvhBK4eXbbg86ioiISKF5KTze8j2FFMjC5X8yxz3PKcl3cNgBlYKOE3NO\nk6ZEROKel5VLu0UjiOxb5guPgyUxpE1JXCpF02lFRBKB5z3UzawBUIvw2h0LnHNf+ZZK/mXm10v4\npvQLnJvSjVo1qgYdR0RExJNCDy41s2pmNhP4AugPDATmmtkMM4v6J6KZdTKzOWb2l5mtMrPxZnZ0\ntK8btBtGdCZpazVG3Xl30FFEREQ88zKrZQDh1UprO+f2c85VAo6LtPX3M9wenB7JcDJwNlAamG5m\nZWNw7UAMnT6HFemvcX2Nx6iSXi7oOCIiIp55udVyPnC2c27hzgbn3AIzux2Y7luyPXDONc373Mxu\nILxnTAPg42hfP9ZCIUf7yR0ok3QcLz7UKug4IiIiReKl8EgCduTTvgNvPShFlUF4nMm6AK4ddd2z\nJrOx0iwePWpSidiTRUREEpuXQmEm0M/MDtzZYGYHAX2BGX4FKwgLb4n7LPCxc25BLK8dC9uyc3hq\nbgcy1p/Jwy0vCDpO4LQ7rYhI/PO6gNgEYKmZrSDc21AD+Ba41sdsBfE8cCxw6r5ObN++Penp6bu0\nZWZmkpmZGaVoRdfuxeFsT1/AS6cOJymppE8nLen//SIi0ZWVlUVWVtYubRs3bvT9Oub1r0gzOwf4\nD+FPhAXOuff8DFaA6w8ELgJOd84t38t59YG5c+fOpX79+FlU9Zc/1nPEs/+hRk4TlvYZHXScwKXc\ndwx1Uy9mzhO9go4iIlJizJs3jwYNGgA0cM7N8+N7FrrHw8zKOOe2O+feBd71I4SHDAOBi4FGeys6\n4lmzvp1wpbbx5s29g44iIiLiGy+3Wjaa2Rzgg8jXp865rX6G2hszex7IBJoDm81s/525nHPbYpUj\nml6a8ikL0wZxRbkB1D/qwH2/QEREJE54GVx6NjAVOIXwWI/1ZvaxmT0euf0Sbe0IrxnyAfB7nq8W\nMbh21G3ZtoO7321LuQ0nMOqeW4OOIyIi4isve7V8THi9jCfNrBRwAtAW6Ah0AqI659M5F8SU3Zhp\n0bcf2yp8z8hGX2j6rIiIJBxPe7WY2VFA4zxfqcAk4H2fcpVIs79fxqRNXakbupNrzoqfgbCxot1p\nRUTin5fBpb8C5YBZhG93PA1847TIQpFdMfgukkpVYtKDjwUdpdgxTacVEUkIXm5brCFceFSLfFUF\nEnaflFjpMPRNVmZM5N5a/TiwcoWg44iIiERFoQsP51xdoDrQk/Atlh7AGjObbWaP+5yvRJj/80p6\nL2rLARsuoecNlwUdR0REJGo8jfFwzm0AJprZx8BswmtqXE14pksX/+IlvlDIcc6AmyEliRntX9IK\npSIiktC8jPG4lH8GldYmvDnbx8D9hMd8SCG06v8yqzMm8fARb1OrRtWg44iIiESVlx6PQcCHwGDg\nA+fct/5GKjnenfsjo/5sz3923EL3ay8MOo6IiEjUeVnHo1o0gpQ027JzuHTkdZRKOoD3O/cJOk5c\n0HRaEZH45+VWSy5wgHNu9W7tlYHVzjmtelUAzZ7qweaKXzDolI+pvl/5oOMUf05jX0REEoGX6bR7\n+gQoA2QXIUuJ8cKk2czM7cZpPMQtFzQMOo6IiEjMFLjHw8zuijx0wM1mtinP4WTgDOAHH7MlpO9+\nWcUds1pQMbsh0556JOg4IiIiMVWYWy3tI/8a4Y3acvMcywaWRtplD7Zl53B6v5ZQJpdZd46lXGrp\noCOJiIjEVIELD+fcYQBm9j5wmXNufdRSJahG3TuzIf0j+jeYSd0jDgg6joiISMx5Wbn0TOfcejNL\nMbNjIjvUyj50HDaeOaV70Sz1Ke5sfkbQcURERAJR6MLDzMqa2SvAFuB7oEakfYCZPehzvoTw7twf\n6bX4Bg7ceBkTO94XdBwREZHAeJnV8hRQh/DKpdvytL8HXOVDpoSyfPVGLhp1KaW3V+fzTkO1JHpR\naANkEZG45+U2ySXAVc65z8ws7yfB98AR/sRKDFu27aD+ky3YXvZX3r7kUw6uWjHoSHFMBZuISCLw\nUnhUBVbn054GWlpyp1DIccIjd7K24kx6153GhSfXCjqSiIhI4LzcavkSaJbn+c5i42bg0yInShDN\ne/ZhYdogbqgyiPsuOyvoOCIiIsWClx6Ph4ApZnZs5PV3m1ltoCHQyM9w8arD0DeZtL0D/xd6iKF3\n3Rh0HBERkWLDy3Taj4G6hIuOb4FzgVVAQ+fcXH/jxZ+h0+fw9JJrOeSvFszq+ljQcURERIoVT2tw\nOOeWAG18zhL3xs/+jptmXkCFbfX4uvtQSiV7uZMlIiKSuDwVHmaWBBwJVGO3XhPn3Ic+5Io7M776\niSsmnEOZHYfwbadJ7FexbNCREo7T2GURkbhX6MLDzE4BRgOH8u85jo7whnElyheLfuW8UWeT7Cry\n5d3TOXT/jKAjJRzTdFoRkYTgpcfjRf6Z2fIHJXwK7cLlf3LaoHMg2fHRTe9Ru2a1oCOJiIgUW14K\nj6OAK5xzP/kdJt788sd6GvQ9j5yU9Uxt8REn1zok6EgiIiLFmpfRj58THt9Roi1asYbavZqwrcwy\nXrtwOuc0OCroSCIiIsWelx6PAUAfM6tOeDrtjrwHnXPz/QhWnH2/dDUn9D+b7SkrGdvsfa44/fig\nI4mIiMQFL4XHG5F/h+Rpc4QHmib84NKvl/zByc81IafUet669AOan3Js0JFERETihpfC4zDfU8SJ\nLxb9yqkvnUVu8hYmt5jFeSccHXSkEsVpd1oRkbhX6MLDObcsGkGKu2lfLubCMefhzDHj2g9pXOfw\noCOVMJpOKyKSCLS0ZgG8Mu1zLhh3Kkmhsnx0o4oOERERr1R47MOjoyZx84dnUX7b0Sy4/2MaHlsj\n6EgiIiJxy9OS6SVF6/5DGLb2Fg7YchHfdR+tZdBFRESKSIVHPkIhR5PHuvMBj3LstnZ81WsgKaUT\nerKOiIhITHguPMwshfw3iVte1FBBWrNxC/W7tWZF+ljOSXqCqT06kZSkgY0iIiJ+8LJJ3FGE1/D4\nv90PEefreHy5+DcavXAxW8ot5IEab9Cr9WVBR5I8tDutiEj889LjMQzIAS4kgTaJGzp9Dje/dwlW\nqhRjzp3NVY3qBh1JdqFeJxGRROCl8KgLNHDO/eB3mKC0eW44g/9oR/kddZl953iOP7x60JFEREQS\nkpfCYwFQxe8gQdiwaRsNu9/ND2kvcdT21sx57HkyyqcGHUtERCRheVnHoyPQy8wam1llM6uY98vv\ngLszs9PNbKKZ/WZmITNr7uX7fDj/Fw58+FR+KDOc6ysNZvHTQ1R0iIiIRJmXHo/3Iv/O2K09VoNL\n04CvCQ9wfWMf5+ar2+jJdPv2WpKTMxh55idcc1Z9XwOKiIhI/rwUHmf6nqIQnHNTgakAZlaoEYd/\nbd7OWU88xNwyz1Bt24V81mEEhx1QKSo5RURE5N+8bBI3KxpBom3290s57dWb2VrhO5qX6cMbve+h\nVLJWjI8nmk4rIhL/PH3yRsZZjDSzT8zsoEjbdWZ2mr/x/HPXB9eQm7SF0Wd9zoQH71XREWdM02lF\nRGIuFPL/D75Cf/qa2eXANGArUB8oEzmUDjzkXzR/1dhxHiu6zCWzcb2go4iIiBR7C5f/SdO+HX3/\nvl7GeHQB2jnnRphZyzztsyPHiqWay9Zz8/WZu7RlZmaSmZm5h1eIiIiUHFlZWWRlZQGwYMUfLPl7\nPmwP+X4dL4XHMcCH+bRvBDKKFid6+vbtS/36mr0iIiKSn8zMTE5pfD5n9b6TpQ3e5oANl/DsWbdx\nVfNzfb2Ol4EOK4Ej82k/Dfi5aHH2zczSzKyOme1c0/zwyPNDon1tERGRRPXwq29zZN/jWJryDm2r\njuDXPm9y5EGVfb+Olx6Pl4F+ZnYj4XU7DjSzhkBvoLuf4fbgBOD9yLUd0CfSPhy4MQbXFxERSRjz\nf15Js4F38Wv661TdcQGTbnqJE485OGrX81J4PEW4p2QGUI7wbZftQG/n3EAfs+UrMp1XU1JERESK\nIBRytB7wCq+ufABKl+aO6qPp93BLkpKiO4vQyzoeDnjCzJ4mfMulPLDAObfJ73AieYX/ryciIkU1\n8bMFtBpzGxsrzeLIHTcw+Z7eHHWw/7dV8uOlxwMA51y2mf0N/KWiQ6LOaR0PEZGiWrluExc/8xhz\nkp6hdHJNetZ+lw5XnB3TDIUuPMysFNAVuItwbwdmtgkYAHRzzu3wNaGIiIgUSSjkuH/oOPr90J5Q\nmbU0KdWVcd3uD2RzVC89HgOAy4AOwKeRtobAo0Bl4FZfkomIiEiRjZ/9HTe/3p51ld6jes7FvH79\ns5x2XM3A8ngpPK4GWjrnpuRpm29my4ExqPAQEREJ3KIVa7h84CN8nzqI0slH0PXId3j0mmZBx/JU\neGwHlubTvhTILkoYERERKZpNW7O5rv/zTNjQDZfkuLhsb0Y+eDvly6YEHQ3wVngMBB42s9bOue0A\nZlYG6Bw5JiIiIjEWCjnuG/I6Axc+RE75X6gdasvrt3ejVo2qQUfbhZfCox7QBPjVzL6JtNUBUoAZ\nZvbmzhOdc5cVPaJImEPTaUVE8tNvwiw6f9CBzRlzqBZqxovnvsWlpx4XdKx8eSk8NgBv7Na2wocs\nInuh6bQiIrsb99F8bn+zM6sz3iGNE+lb533uuaRx0LH2yssCYq2jEUREREQKZtqXi7klqyvLK7xG\nqaTDueeg1+jz8JVRX3XUD54XEDOzqoR3qnXAYufcn76lEhERkX/5dMFybhjancXlhpGcdADXZQzi\nxQdvoFxq6aCjFZiXBcTSCK/l0Yp/9kzJNbMRwJ3OuS0+5hMRESnxPl2wnJuH92BBmVewpAwuKdeb\noQ+0C2QBsKLy0uPxDNAIuAiYHWk7DehPeKdYreMhIiLig9nfL6PN8B4sTB2CWToXlHmcIffcRvX9\nygcdzTMvhcflwBXOuQ/ytE02s63AWFR4iIiIFMmMr37ijqxe/JA6DEtKp2nqE7xyz61xXXDs5KXw\nKAesyqd9deSYSJRoOq2IJLbxs7+j/Zs9WFZ+DElWNaEKjp28FB6fAt3MrJVzbhuAmZUlvHHcp3t9\npYhHpum0IpLABk/9jIenPcXKjAkkJ9XgqooDeP6+1uxXsWzQ0XznpfC4G5jKPwuIOaAusA04z8ds\nIiIiCSsnN0T3rMn0+7IXf1X6iJSkY7ip8lD6d7wmrmapFJaXdTy+M7OjgGuB/xBe2WkMMMo5t9Xn\nfCIiIgnlr83buXdoFiN/fprt6Qsob6fwYM3xPHZtc0olJ+37G8Q5T+t4RAqMl33OIiIikrAWrVjD\nHcMGMfPvgYTSVlLNXUjXBoNo1/TUuFj4yy+eFxATERGRfZs85wceeKMfC0oPBxy1uJ7eF9xD05P+\nE3S0QKjwEBER8VlObognx06j/5x+rM2YRpKrTpOUh3iudTuOOaRK0PECpcJD4oZzmk4rIsXb72v/\n5u4hw5mwcgA7Ki6mHA1oW3UEve5vQcW0MkHHKxZUeEicKDn3P0Uk/kz8bAGdxz/Hd6VGQKmtHOIu\np/OJQ2lzfsMSNX6jIIqySVwDoBbh6bQLnXPzfEslIiJSzG3ZtoOuoycyeP5zbKj0Pknsz+ml2tPn\n6ls48ZiDg45XbHnZJK4a4emzjYENhP8UTTez94GW2qVWREQS2ecLV/DAmJeZvXUwobQ/qMCp3FF9\nND1aXU75silBxyv2vPR4DAAqArWdcwsBzOxYYDjhjeIy/YsnIiISvOwduTw1bjovzHmRlRXfgZxy\n1LbreOTMtrQ4o07Q8eKKl8LjfODsnUUHgHNugZndDkz3LZmIiEjAvlj0Kx3HDOHDTa+QW345qdTh\n6oznebrV1RxYuULQ8eKSl8IjCdiRT/uOyDEREZG4tS07hydem8LL815mVcVJkFOWY8ikY8M2XH/2\niRosWkReCo+ZQD8zy3TO/Q5gZgcBfYEZfoYTyctpd1oRiaIZX/3Ew28O4fPtwwil/UFZ6nF1xvP0\nvDaTg6tWDDpewvBSeNwBTACWmtkKwrNaagDfEt6/RSQK9BeGiPhvzcYtdBn1Bq8tHsKGSh9ATjr/\nTbqWTo1uIrNxvaDjJSQvm8StAOqb2Tn8s0ncAufce36HExER8Vso5Bg87TP6vj+UH0qNgTJ/k0Fj\nbt1/JI9ffVlCbkVfnHiZTtsKeM059y7wbp72FMLTaUf4mE9ERMQXXy7+jYdfH8nMtcPITv+BZFeD\nM1La0+2S62lc5/Cg45UYXm61DAWmAqt3a68QOabCQ0REioV1f23l0TETGL1wGGsrvgu5KdS0y2l7\nbH/uv6xJidiGvrjxUngY5DvK72BgY9HiiIiIFE0o5Hhx8mz6zxrOotJjocxfVHCn0mq/QTyWeSU1\nqqUHHbFEK3DhYWZfES44HDDDzHLyHE4GDiPcEyIiIhJzM79eQve3XmX2plfJqfAzye5QTit9F12a\nXcd5JxwddDyJKEyPx1uRf+sC04BNeY5lA0uBN/yJJSIism9Lfl9Hl6zXeGfFq2yq9ClkV+Aou5Lb\n6w7h9gtP162UYqjAhYdzrhuAmS0lPLh0W7RCieRH63iICMCGTdt48vXJjJz/Kn9UmAQWoqqdyw3V\nR9O15cVUSS8XdETZCy/TaYdHI4jI3pjW8RAp0XJyQ7ww6WOe/3gki5Jfx6VuoBwNuCStF4+1yOS4\nw/YPOqIUkJfBpSIiIjHxxsff8vTUUXy5PYvc8sspFarJqWXvoNMF19D0pP8EHU88UOEhIiLFyqcL\nltN9fBYfrB3FtvRvsZzKHJvUgttPuIa2F/yf9kqJcyo8REQkcItWrKHb6+OYtHw0f1X6CHaUpYZd\nzPVHPEmHy8+lfNmUoCOKTzwXHpGVSg8DljjncvZ1vt/M7HbgfqA68A1wp3Pui1jnEBERb1au28Rj\nYycyblEWqytMBXNUsXNoW3UEj1x1ibadT1BelkwvBwwAro80HQ38bGYDgN+cc0/5mG9PGa4C+gC3\nAHOA9sA0MzvaObcm2tcXERFvNm3Nptcb0xnx1WiWpU6AlC2U5xSuqNCXrldcqUGiJYCXHo8eQB2g\nMbsuGPYe8CgQ9cKDcKExaOe+MGbWDmgG3Aj0isH1JQDOaTqtSDzKyQ0x8O0PGfRJFouSxuHKrqOM\nO5ZzUjvzUPOW2ielhPFSeFwCXOWc+8zM8n4SfA8c4U+sPTOz0kAD4Mmdbc45Z2bvAQ2jfX0JigaT\nicSTUMgx+v15PDtjNF/teI1Q+d9IDh3KyWVu4b6zMrns1P9qkGgJ5aXwqMq/N4gDSCP/PVz8VoXw\nEu2rdmtfBRwTg+uLiMgeTPliET0nZfHJ36PZUfFHLLcqxyVfxe0nZtLm/IYqNsRT4fEl4dsaAyLP\ndxYbNwOf+hHKoz1tXiciIlH0+cIVPPbWa8xcncXWjHmQXZEj7TJurDWQ9pecRWqKJlDKP7z8v+Eh\nYIqZHRt5/d1mVpvwbY5GfobbgzVALrD7CKRq/LsX5H/at29PevquOxJmZmaSmZnpe0ARkUS3c/rr\n5OVZbKz0IeSU4SC7kKtrdOahK5uSUT416IhSSFlZWWRlZe3StnGj/5vOe1ky/WMzqws8CHwLnAvM\nAxo65771OV9+199hZnOBJsBEADOzyPP+e3pd3759qV+/frTjiYgkrJXrNtF97ATGLRrNnxWmgzkq\n29ncXGUYD195ibabj3P5/TE+b948GjRo4Ot1PPV/OeeWAG18TVI4zwDDIwXIzum05YBhAWYSEUk4\nm7Zm02PcVF79ajQryk2E0lup4P6PKys+S9crrqR2zWpBR5Q442Udj3eBkcCbzrm//Y+0b865sWZW\nBehO+JbL18B5zrk/g8gjsaHdaUVi45/pr6NZlDwOl7qeVP7LeamP0OWSlpx2XM2gI0oc89LjsZDw\nWh0vmNk7hIuQKc65Hb4m2wfn3PPA87G8pgRHu9OKRFco5Hjtw6955t3RzMvOIlT+N0qFatIw9Vbu\nPzuTS089LuiIkiC8jPG4y8zuAc4GrgZeBXLMbBww2jk3y+eMIiISJR988zOPTxjNRxtGk52+EMut\nQu2kFtxx4jWa/ipR4XWMRwiYDkyPrBp6EdCZ8JTaZP/iiYiI3xatWEPXsWOZvGIUf1f6BLLTOMwu\n5Yaj+3D/pWdTLrV00BElgRVpcrWZVQdaAtcCxwPapE1EpBha99dWur02kTELRv5vQ7aqdh7XVx/N\nwy2aU61SWtARpYTwMri0InA54dssjYGfgdFAS+fcT76mExERz7J35DLg7VkM+nQkP5YeB2X+pjyn\ncGXFZ+nWogW1alQNOqKUQF56PFYB64GxwEPail5EpHiZ8Mn39Jj8Kl9uH0Vu+V8pFTqCxmXuo0vz\na2hS78ig40kJ56XwuBh4LzLOQySGNJ1WZE++X7qah8dmMfWPEWzNmIfl7MexSVdx10nXcfN5p2iQ\nqBQbXma1TI9GEJG90y9Nkd39tXk73ca8zcjvhrO64hRwSRzAhdxd82E6t2hK+bIpQUcU+ZcCFR5m\nNg9o4pxbb2ZfsZc/PZ1zWpdcRCRKQiHH8Pe+oM+MYSxIGoNLXU8aJ3FlhX480bIlRx1cOeiIIntV\n0B6PCcD2PI/V5y0iEkPzf15JpzGvMmPNMLanLyApdBAnp7Sl8/mtuPDkWkHHEymwAhUezrlueR4/\nGrU0IiLyP5u2ZtN9zDsMnz+E1RWnQqgUNexS2tXqy32XNiGltJZNkvjjZTrtz8CJzrm1u7VnAPOc\nc4f7FU4OgMCWAAAgAElEQVREpCSa8Mn3PPbOEOblvoor9ydpnERm+kB6XNOSQ/fPCDqeSJF4mdVS\nk/xXJy0DHFykNCIiJdTKdZvo+OoY3lg6mM0Zn2OhKtRLvo4uZ9+ofVIkoRS48DCz5nmenmdmG/M8\nTwaaAL/4FUxEJNHtHCja892XWZQyBkpvpirnc+sh4+iaeZFmpUhCKkyPx1uRfx0wfLdjO4ClwH0+\nZBLJl3Ma0yyJYfnqjXR4dRQTfh3Etoz5JFODxin381SLGzm51iFBxxOJqgIXHs65JAAz+4XwGI81\nUUsl8i9ax0PiWyjkeHXGlzw5fRCLU7Kg1HYOoDm3HdOTDpefo4GiUmJ4WUDssGgEERFJRGs2buH+\n4Vm8/ssLbMmYSzKH0KRMJ3pl3kj9ow4MOp5IzHnandbM0oBGQA1gl5uQzrn+PuQSEYlrU75YROc3\nX+RrhuHKbKQaTbnviLfpctUF6t2QEs3LdNp6wGSgHJAGrAOqAFuA1YAKDxEpkbJ35PLYmMk8P3cA\n6yq9i7kqnFSqLb0ua8sZx6uzWAS89Xj0Bd4G2gEbgVMIDy4dCfTzL5qISHz45Y/13DVsCFPXPkdO\nhV9IsxNpW3UET7W6kozyqUHHEylWvBQedYG2zrmQmeUCZZxzP5tZB8KzXd70NaGISDE15YtFPDCu\nH9+XGg7JOzjMXcXDDcfQ+tyTgo4mUmx5KTx28M9eLasJj/NYSLj3o4ZPuUT+xWmLICkGQiFHrzfe\n45lPnuXPjMkksT9npnSkf6u2HHfY/kHHEyn2vBQeXwEnAIuBWUB3M6sCXAd862M2kf8xTaeVgG3a\nmk37IVm8uqQ329O/oyx1ubnKMPrc0JKKaWWCjicSN7wUHg8BFSKPOwMjgBeAH4EbfcolIlIsLFu1\ngXaDBzF9Y39Cab9TzV1IpzoDuKt5I5KSVBCLFJaXdTy+zPN4NXC+r4lERIqBLxf/RrthzzDXXoKk\nbI5xreh9/r3agl6kiDyt4yEikqimfLGIe8Y+zeKyI8CV4/9K3cULN9zJ8YdXDzqaSEIoUOFhZl9B\nwUb2OefqFymRiEgAXpv1NfdNeILfKr5BklWnWdknefHOWzi4asWgo4kklIL2eLy171NEROLP0Olz\n6DTlcVZlvE2ppMO5LmMQAx9opQGjIlFSoMLDOdct2kFE9kXTacVPgyZ/Qufp3VhbaTopScfQtuoI\nnu2USWqK7kCLRJPXvVoygCuAI4CnnXPrzKw+sMo595ufAUUAcJo9IP4YOn0OD0x5hLUZ0yiTdBx3\nHTiGpx+6QvuniMSIl71ajgfeI7xgWE3gZcL7tVxGeAGxVj7mExHxxaiZ87j37a6szniHFKtF+4PH\n0qvL5ZRKTgo6mkiJ4qXH4xlgmHOug5n9nad9MjDan1giIv6Y8sUi2mR14bf0cZROOprb9h9F34eu\nUg+HSEC8FB4nAm3zaf8N0HwzESkWvlj0K9cN7saickNJTj6ImyoPZWCnazWGQyRgXn4CtwP5zS87\nGvizaHFERIrmlz/Wc9XAHnxhA7Dk8lya1odhHdpplopIMeGl8JgIPGJmLSLPnZnVAHoCb/iWTESk\nEDZtzeb6AS8yfl03XPJ2GpXqyMi779U6HCLFjJfC4z5gHOGdacsS3iiuOvAp4b1bRKJE02nl30Ih\nR6cRb9H3247sKL+EY3JvZEyb7tQ94oCgo4lIPrzs1bIROMfMTgXqAOWBec659/wOJ7KTdqeV/Iz7\naD43v3EXGyvNonLuuQw65w0uP+2/QccSkb0oVOFhZqWBqUA759xsYHZUUomI7MWS39dxaf9H+LbM\nC6QkH0W3oyfzSOYFQccSkQIoVOHhnNsRWcdDRCTmsnfk0nrgYLJWdcYlZXNhai+yHryT8mVTgo4m\nIgXkZYzHSOAm4EGfs4iI7NHIGXNp+05btmTM5cicG3jjth7aMVYkDnkpPEoBN5rZOcCXwOa8B51z\n9/oRTEQE4Nc//+LCZx7mm5SBpFptXjxxNm2b/l/QsUTEIy+Fx3HAvMjjo3c7pmkHIuKLUMjRcdib\nPPPDXYRKb6BZak/Gdrqbcqmlg44mIkXgZVbLmdEIUhBm9hDQDKgLbHfO7RdUFok951TXlhTzf17J\n+QNu44+M8VTf0Zw3bxhAw2NrBB1LRHwQb7sjlQbGAi8EHURiTdNpS4JQyHHLcyOo+/KxrEyZTfuD\nx/JH3wkqOkQSSFxtWuCc6wZgZtcHnUVE/PXFol9p9mJb/syYTM3sa5h6z7Mcc0iVoGOJiM/iqvAQ\nkcQTCjnueGk0Lyy7naTSaXQ5fCKPXXdR0LFEJEpUeIhIYJb8vo4ze9/KivSx1My+hpn3D+CwAyoF\nHUtEoijwwsPMegAd93KKA2o55xYX5Trt27cnPT19l7bMzEwyMzOL8m1FxKNe496j0+c34Mps5q4D\nx9Cv61VBRxIp0bKyssjKytqlbePGjb5fx1PhYWbXAe2Aw4CGzrllZnYP8ItzbkIhv11vYOg+zvnZ\nQ8xd9O3bl/r16xf124hIEW3amk3jxx9ibkofKu1owrSbh3HiMQcHHUukxMvvj/F58+bRoEEDX69T\n6MLDzG4FugPPEt6NNjlyaANwD1CowsM5txZYW9gcIhJ/Ppz/C02HtmRzha+4uMwzjOtzN6WS421y\nnYgUhZef+DuBNs65J4DcPO1fAlHdFtLMDjGzOsChQLKZ1Yl8pUXzulI8OK1PF9c6DhtPo6x6bE9e\nzbDTZ/PWg+1VdIiUQF5utRwGfJVP+3Yg2gVAd6BVnuc7V1A9E/gwyteWAJnW8Yhbm7Zmc1r3B/gm\ntT8Hbr+MTzq+wqH7ZwQdS0QC4uXPjV8Irxy6u/OBhUWLs3fOudbOueR8vlR0iBRDXy/5gwM7nck3\npV7kinIDWNF7nIoOkRLOS4/HM8BzZpZKeDnJk8wsE+gE3OxnOBGJXy9N+ZRb378cUoyXT53Fzeef\nEnQkESkGvOzVMtjMtgKPA+WA0cBvwN3OuTE+5xOROHTtsy8xat0dVMg+iY/vGqft60XkfzxNp3XO\njQJGmVk5oLxzbrW/sUQkHm3ZtoMTHrmThWmDqL39Vj7r8Szly6YEHUtEipFCj/Ews5lmlgHgnNuy\ns+gws4pmNtPvgCISH5av3sghDzZlYeoQWmW8zHc9n1fRISL/4qXHozGQ32+TVOD0IqUR2QtNpy2+\nZn+/jCavNCO77G/0qTeNey89M+hIIlJMFbjwMLPj8zw91szy3rRNJjyr5Te/gonsStNpi6sR731J\n6+kXkZRUlomXfMKFJ9cKOpKIFGOF6fH4mvC+KQ7I75bKVsKLi4lICdH51Yk8uaglaTuO5/O7J1K7\nZrWgI4lIMVeYwuMwwn92/gycBPyZ51g2sNo5l5vfC0Uk8dzy3AheXn0jB265mG+7jWS/imWDjiQi\ncaDAhYdzblnkodY4FinhLu/Vnze33s3RW27i256DSCmdvO8XiYjgcTotgJkdC9Rgt4GmzrmJRQ0l\nIsVTKORo8lh3PuBRTsi+n8979iIpSeNvRKTgvOxOezgwnvCGcI5/Rv3tnHKgP31EElBObogTu9zL\n16n9ODf5SaY89qCKDhEpNC+3TfoR3q9lf2ALUBs4g/DutI19SybyL5pOG5Sc3BDHd7qNr8v0p2X5\n55nWpZOKDhHxxMutlobAWc65P80sBISccx+bWSegP1DP14QiaHfaIIVCjnoP3cnCci/RuvIrDLmz\nddCRRCSOeenxSAY2RR6vAQ6MPF4GHONHKBEpHkIhR/3O9/BduedpVeklFR0iUmReejy+A44nPK32\nc6CDmWUDt0TaRCQBhEKOE7vcxzep/bm64osMv1ubT4tI0XkpPB4H0iKPHwHeAT4C1gJX+ZRLRAIU\nCjkaPvIg88r05cq0gYxq3zboSCKSIApdeDjnpuV5/BPwHzPbD1jvnNPoP5EEcGGPp5lTuheXpPZl\n7P23Bx1HRBJIocZ4mFkpM8sxs+Pytjvn1qnoEEkMtzw3gik5HTkt1IXxHe8JOo6IJJhC9Xg453LM\nbDlaq0MCoNo2+h4fM5WXV93E0VtvYlbP7kHHEZEE5GVWyxPAk5HbKyIxoum00Tb83S94+Nsr2P/v\nC/jmiRe1ToeIRIWXwaV3AEcCv5vZMmBz3oPOufp+BBOR2Hl37o+0ntGU8tuO57tuY0hN8bybgojI\nXnn57fKW7ylEJDC//LGeZqObUtpVYd79b1MlvVzQkUQkgXmZ1dItGkFEJPa2ZedwQs+ryEldx8zM\nLzjq4MpBRxKRBKf+VJES7NRHO7AufSa960yncZ3Dg44jIiWACg+REqrNc8OZV6Yvl5frz32XnRV0\nHBEpIbzMahEJhNPutL4ZPPUzBq+8haM33cTY++4IOo6IlCDq8ZC4oN1p/fPl4t9oO/NSKmSfwBdP\nPKdpsyISU4Xu8TCzxlHIISIxsC07hybPZ4JL5pO736RiWpmgI4lICePlVss0M1tiZl3M7BDfE4lI\n1Jz/5OP8lT6bfmdkcdxh+wcdR0RKIC+Fx0HAQOAK4Bczm2ZmLcwsxd9oIuKnZ9/6gFnuMc5KfpQ7\nLjo96DgiUkIVuvBwzq1xzvV1ztUFTgIWA88Df5hZfzOr43dIESmahcv/5L5PriFjYyOmdHoo6Dgi\nUoIVaVaLc24e0INwD0gacCMw18w+MrPaPuQTkSIKhRyN+t6AS8pmxu0jSSmtPR5FJDieCg8zK21m\nV5jZZGAZcB7hPVz2J7yPy6/A676lFBHPLu3Vlz8zJtP1+OHUP+rAoOOISAlX6Om0ZjYAyIw8HQl0\ncM59l+eUzWbWEVha9Hgi/9A6HoU39sNvmLjlQRrk3kfXq5sGHUdExNM6HscCdwJvOOey93DO78CZ\nnlOJ/IvWmiisLdt2cMP41pSxY/jgiSeDjiMiAnjbJK5JAc7JAWZ5SiQivri0d2+2VvyG4ad/Tvmy\nmnQmIsWD55VLzexYoAawy28059zEooYSkaKZ+NkCpm9/lJPdA7Q6+4Sg44iI/I+XMR6HA+OB/wKO\nf/rAd96A15B5kQBl78jl6tduJCXpMCZ37Rp0HBGRXXiZ1dIP+IXwDJYtQG3gDOBLoLFvyUTEkyv7\n9GNz+hz6NxnCfhXLBh1HRGQXXm61NATOcs79aWYhIOSc+9jMOgH9gXq+JhSRAnt37o9M3NSZerl3\n07bp/wUdR0TkX7z0eCQDmyKP1wA7FwZYBhzjR6g9MbNDzWywmf1sZlvM7Ecze9TMSkfzulJcaDrt\n3oRCjhav3kKpbQcytcPjQccREcmXlx6P74DjgZ+Bz4EOZpYN3BJpi6b/EB5T0gZYAhwHDAbKAR2i\nfG0JkGk67T61f2UsGyp9wJO1plGtUlrQcURE8uWl8Hic8PLoAI8A7wAfAWuBq3zKlS/n3DRgWp6m\npWbWG2iHCg8pwdZs3MJzix+gek5zOrU4N+g4IiJ75GUdj2l5Hv8E/MfM9gPWO+eC6AvPANYFcF2R\nYuPyvj3JTV3FqMufCTqKiMheFXqMh5mlRwqN/3HOrQMqmVlF35IVLMuRhPeIeTGW1xUpTmZ/v4wP\nc3pxCvdyVt0jgo4jIrJXXgaXjgFa5tPeInKs0Mysh5mF9vKVa2ZH7/aag4ApwGvOuSFeriuSCK56\n5X6Ssivxxj3a7l5Eij8vYzxOBu7Np/0D4AmPOXoDQ/dxzv8GrprZgcBM4GPnXNuCXKB9+/akp6fv\n0paZmUlmZuYeXiFS/D0z/n1+Sx9H26ojOLByhaDjiEgcy8rKIisra5e2jRs3+n4dK+ywDDPbDJzi\nnPt2t/b/Ap8758r5mC+/6x9EuOj4ArhuX+NKzKw+MHfu3LnUr18/mtEkiva7pwkVkquxrE/Wvk8u\nIbZl55DxYD1Ku/Ks7z2bUsleOjBFRPZs3rx5NGjQAKCBc26eH9/Ty2+qOYSnzu6uHTC3aHH2zswO\nINyzspzwLJZqZra/me0fzetKcaDptLtrPeBltqd/x8Cm/VV0iEjc8HKrpQvwnpnVAWZE2poAJwLR\nnsd3LnB45GtFpM0IryylPWKkxFj311bGrnqMmrnXcP05JwYdR0SkwAr9Z5JzbjbhZdN/JTyg9CLg\nJ+B459xH/sb717WHO+eSd/tKcs6p6JAS5cbnXyRUdjWDr3006CgiIoXipccD59zXwNU+ZxGRAli9\nfjMT1z/FUaHraVLvyKDjiIgUiqcbw2Z2hJk9bmajzaxapO0CM6vtbzwR2d0Nzz+HK7OOwdd1CTqK\niEiheVlArBHwLeFptZcD5SOH6gDd/IsmIrv7fe3fTP27F7W23cQZxx8WdBwRkULz0uPxFNDFOXcO\nkJ2nfSbhsR8iURHMivzFS6vn+uFK/83Q1p2DjiIi4omXwuO/wPh82lcDlYsWRyR/2p0Wlq3awMyt\nfTh+R1tOrnVI0HFERDzxUnhsAA7Ip70e8FvR4ojInrR6oS8ueRvD23QKOoqIiGde92rpaWbVCa+f\nkWRmpxJe9nyEn+FEJGzJ7+v4MLsvDdzt1D0iv7pfRCQ+eCk8HgJ+ILyAV3lgAfAh8AnwuH/RRGSn\n24YMguRshrfpEHQUEZEiKfQ6Hs65bKCNmT0GHEe4+PjKOfej3+FEBDZtzea9vwbwH1pRu2a1oOOI\niBSJpwXEAJxzywnvmSIiUXT/0NcIpf1Bz/PuCTqKiEiRFajwMLNnCvoNnXP3eo8jsjclbzptKOQY\n8WNfqnA+zU85Nug4IiJFVtAej3q7PW9AeFO2RZHnRwO5RHl3WinJSuZ02v4TZ7E14yserjUt6Cgi\nIr4oUOHhnDtz52Mzuxf4G7jeObc+0lYJGApEdZM4kZKm56y+lLHadLzinKCjiIj4wsuslvuATjuL\nDoDI4y6RYyLig3fn/sjK9LdpWbM9SUkls8dHRBKPl8GlFYGq+bRXBSoULY6I7HTf6/0wV4Vnb7wm\n6CgiIr7x0uMxHhhqZpeZ2cFmdpCZXQ68ArzpbzyRkmnJ7+v4NnkojcrdTkb51KDjiIj4xkuPRzvC\nq5SOBkpH2nIIFx4P+JRLpES7Y+jLkJTLwOtvDTqKiIivvCwgtgW4zcweAI4gPN3gJ+fcZr/DiZRE\n2TtyeXfDQI7iGi0YJiIJpygLiG0G5vuYRWSvXAlZx6PXG++SW/5XOv9fu6CjiIj4znPhIRJLVoLW\n8Rg0ZwhlqM11TU4IOoqIiO+8DC4VkShZ8vs6fk2bwHnVbtQUWhFJSCo8RIqRB0eNBgvxZMtrg44i\nIhIVutUiUoxM+n0o1WmmQaUikrDU4yFSTIz98Bu2Zsyjdb3WQUcREYka9XiIFBM9pgzFQtXo0qJp\n0FFERKJGPR4SNxJ5Ou2mrdl840ZRP/layqWW3vcLRETilHo8JD5YYs/w6D7mHVzZNXQ+W7dZRCSx\nqfAQKQaGzx9KGidy6anHBR1FRCSqdKtFJGDzf17J6opTuOhg9XaISOJT4SESsE5jXoVQKZ68umXQ\nUUREok63WkQC9sGfYznILuKwAyoFHUVEJOrU4yESoE8XLGdLxpdcesxlQUcREYkJFR4SRxJvOu3T\n74yHnBQ6XtYs6CgiIjGhWy0SFxJ1d9qZv79JVTubg6tWDDqKiEhMqMdDJCDf/bKKjRkf0fQw3WYR\nkZJDhYdIQHpNmAjO6HhJ86CjiIjEjG61iARkyrI3yXBnUKtG1aCjiIjEjHo8RAKwbNUG1lSYQZOD\ndJtFREoWFR4iAeg5fhIk7+CBCy8JOoqISEzpVovEjUTanXbi4vGkuRM5udYhQUcREYkp9XhIXEik\n6bRrNm7ht7JTOK2ybrOISMkTd4WHmU0ws2VmttXMfjezEWZ2QNC5JHaysrKCjlAkfd6aDilbaH++\nCo+d4v09lV3p/ZS9ibvCA5gJXAkcDVwGHAG8Hmgiial4/6U29ts3KbOxNuedcHTQUYqNeH9PZVd6\nP2Vv4m6Mh3OuX56nK8zsKWC8mSU753KDyiVSEJu2ZvNL6bc5LeXOoKOIiAQi7gqPvMxsP+AaYLaK\nDokHz0/6EJe6gVsbXRp0FBGRQMTjrRbM7Ckz2wSsAQ4BNCdR4sK4r6aTtLk6V51RN+goIiKBKBY9\nHmbWA+i4l1McUMs5tzjyvBcwGDgU6Aq8Cly4h9emAixcuNCfsBKIHas34tjGvHnz2LhxI/PmzQs6\nkifzl06kutXj66+/CjpKsRLP76n8m97PxJHnszPVr+9pzgW/NoKZVQYq7+O0n51zOfm89iBgBdDQ\nOfd5PsevBkb5ElRERKRkusY5N9qPb1Qsejycc2uBtR5fnhz5t8wejk8jPA5kKbDN4zVERERKolSg\nJuHPUl8Uix6PgjKzE4GTgI+B9cCRQHegKnCcc25HgPFERERkH+JtcOlWwmt3vAf8ALwMfA00VtEh\nIiJS/MVVj4eIiIjEt3jr8RAREZE4psJDREREYiYhCg8zu93MfolsHPdZZBDq3s6/0swWRs7/xswu\niFVW2bfCvJ9mdr2ZhcwsN/JvyMy2xDKv7JmZnW5mE83st8h707wAr2lsZnPNbJuZLTaz62ORVfat\nsO+nmTXK83MZyvOzWi1WmWXPzKyTmc0xs7/MbJWZjTezfW4iVdTP0LgvPMzsKqAP4YXE6gHfANPM\nrMoezm8IjCY8MLUu8BbwlpkdG5vEsjeFfT8jNgLV83wdGu2cUmBphAeA3054IcC9MrOawDvADKAO\n0A8YbGbnRC+iFEKh3s8IBxzFPz+fBzjnVkcnnhTS6cAA4GTgbKA0MN3Myu7pBX58hsb94FIz+wz4\n3Dl3d+S5EV5QrL9zrlc+548Byjnnmudp+xT4yjl3W4xiyx54eD+vB/o65/aLbVIpLDMLAZc45ybu\n5ZyewAXOuePztGUB6c65pjGIKQVUwPezEeEdxSs55/6KWTjxJPIH3mrgDOfcx3s4p8ifoXHd42Fm\npYEGhP86AsCFK6n3gIZ7eFnDyPG8pu3lfIkRj+8nQHkzW2pmy81MvVfx7RT085loDPjazH43s+lm\n9n9BB5I9yiDcQ7VuL+cU+TM0rgsPoArhlUtX7da+inCXXn6qF/J8iR0v7+ci4EagOeEVapOATyJL\n6Uv82dPPZ0Uz29PqxFJ8/QG0BS4nvAbTCuADM9MuicVMpHf5WeBj59yCvZxa5M/QYrFkehQYBb//\n6OV8ia09vj/Ouc+Az/53YrjLbyFwC+FxIhL/LPKvfkbjTGRjz8V5mj4zsyOA9oAGDRcvzwPHAqd6\neG2hPkPjvcdjDZAL7L9bezX+XZHttLKQ50vseHk/dxHZSPArwsvpS/zZ08/nX8657ADyiP/moJ/P\nYsXMBgJNCa8C/sc+Ti/yZ2hcFx6RZdLnAk12tkW6i5oAn+zhZZ/mPT/inEi7BMjj+7kLM0sCjiPc\nxSvxJ7+fz3PRz2ciqYt+PouNSNFxMXCmc255AV5S5M/QRLjV8gww3MzmEq6k2wPlgGEAZjYC+NU5\n91Dk/H7ALDO7F5gEZBIe0Ngmxrklf4V6P83sYcK3Wn4iPDCqA+HptINjnlz+xczSCP91u/N2yeFm\nVgdY55xbYWY9gAOdczu73V8E7ojMbhlC+BfcFYT/GpOAFfb9NLO7gV+A7wnvctoGOJPwB5UEzMye\nJ/wZ2BzYbGY7ezI2Oue2Rc4ZDvzm52do3BcezrmxkSlA3Ql3/3wNnOec+zNyysFATp7zPzWzTOCJ\nyNePwMX7GEwjMVLY9xOoBLxEeGDTesI9Jg2dcz/ELrXsxQnA+4Tv/zrCa7QADCc8KLg6cMjOk51z\nS82sGeEC9C7gV+Am59zuo+glGIV6P4GUyDkHAluA+UAT59yHsQose9WO8Pv4wW7trYERkceHEL4F\nDvjzGRr363iIiIhI/IjrMR4iIiISX1R4iIiISMyo8BAREZGYUeEhIiIiMaPCQ0RERGJGhYeIiIjE\njAoPERERiRkVHiIiIhIzKjxEREQkZlR4iIiISMyo8BAREZGYUeEhEsfM7H0zeyZa5xd3Qfz3mFll\nM1tlZjVied09ZNnlv9/MxphZ+yAzieyLCg+RQjCzoWYWyvO1xsymmNl/85xTxcxeMLNlZrbNzP6I\nnNMwCpEuBR6O4vnyb52Bt5xzy4MOko/HgC5mViHoICJ7osJDpPCmAPsT3gL8LCAHeDvP8TeBOsB1\nwFHARYS3na7sdxDn3Abn3OZonV9SmFk5M+tgZi3MrOdezitLeMvwwT5eu7Rf38s59z2wBLjWr+8p\n4jcVHiKFt90596dzbrVzbj7QEzgk0gWfDpwGdHTOfeicW+Gc+9I519M5986evmGky7y/mfU1s3Vm\nttLMbop8IA4xs7/M7EczOz+f1z0TeVwl0rvyYJ7jDc1su5mdufv5eZ73M7OeZrY28vqueY6XN7NR\nZrbJzH4zs3sKcnvDzM4zs4/MbH2kV+htMzu8oNf1em0L62RmP5vZFjP7yswu31vWiNeAd5z7/3bu\nLVaq6o7j+PeHN0SCQmJNNMRUiDEEjNV6rK1Ek7bG0JpgS7w9kLax1BeNJWKbJqUSlSYWtYEHNURU\nIK3RNNRLvEZPTC9pbaEWlGoq+FCvAaV4wCpE/j7815HdcfbsmSPMQf19kp2zZ8/a679mzTnZ/73W\n2ifuBWZIOrGm3LfI7/9vXdRZ18ZBScvKd70FeLQcb+qzcZJWShoq/TG/JsSDwMUjbZ/Z/ubEw+wT\nkDSevLv8d0S8Bewo22xJh/ZY3VxgC3A6sBS4DbgP+BPwJeBxYKWkse1OjoitwA+ARZJOlXQEsApY\nGhGDDXF3AAPANcBCSV8v790CnAl8G/gmMBM4tYvPcgRwE3AaOSr0AbCmh7gjjf0z8vuYB0wrdayS\nNLPuBEmXAWMjYmM5NA44oab4WcDfG9rQjbnA+8BXgcvLsaY+W0L2wfnAucA5pWyrZ4CBfTmSYrZP\nRYQ3b9663IA7gd3AUNn2AK8Ap1TKXABsBd4F/gjcAMxoqHcQeLryekyp/67KsWNKvIGW825uqWsZ\n8JH1uiMAAASNSURBVAKwGngWOKSufGvccuyvwGJgPHlxvKDy3gQyWbi50+dp8/mOLm2f1hS37HcV\nu/p5gEPL+2e01LscWN2hbS8Bl1T6fVtrHZWya4DlLcd+AWwAXgY2Vn7eVVPHILC2lz4jk5L3gO9U\n3p8I7Gzz/c8gk5bJo/334s1bu80jHma9ewo4mVzHMUCORDwqaTJARKwBjiXvTB8BzgbWSZrbUO/6\n4Z2I2AO8RV7Qho+9WXa/0FDPAuBgYA5waUTs7jZu8XqJcUKp56NphYh4B3ixoT4kTZX0G0mbJG0H\nNgMBVJ8EqYvLCGNPJUcrnijTEUOShsi1NlNq2vllYDIwVdJPyBGSQ4B/1sQ4nEwAhs//LvAk+fvw\n64iYBtwWEdMi4nsd2vqxUZOGPptS2vXMcPmI2Eb7/vgfILIvzA44B492A8w+hXZGxMtlf3MZqt8O\n/BBYCBARu8gL0pPADZKWA4uAlR3qbU0Qos0xaJ4inUImPmOAL5J33520izuGvHipvK5SQ30AD5F3\n/pcBr5X6nidHJZriVmP0Ent8+TmrxKx6v+ac04E/R8R1AJIWAr+NiPdqym8lRxqycRG/K+edC6yX\nNJ0cdWnSboFva58dBDxH9lldf7QzqZTb0kVZs77ziIfZvrGHvBuu8y9yuHy/KvP6q4F7yMdmV0g6\neoTVbSKTg4FK/RPIJ3U6tWEScCJwfUQMRsSL9P5EzybyaaFeYm8kE4zjI2Jzy/ZqzTkTyOkoJI0h\nR0eWdojxD3Lqo9X3yWm1tiMrTWr6bFKlyEtkf3ylcs7Eck6r6cArEfH2SNpitr95xMOsd4dJOqbs\nTwSuIJOKB8oF5D5gBTmVMETeVS8Aft+Hti0mL6ZXkGtMZpW2nN9rRRGxQ9LdwBJJ28g76GvJ9QOd\n7ry3kdNE8yS9ARwP/LLhnE8cu5yzBLhF0kFkInAk8DVge0SsanPaJrK/IEesVkfEhjblhj0GLJZ0\nZERsB5B0EnB4ROyWdBg5Bderjn0WETsl3QH8StLbZH9cT/ZHq5nk9J/ZAcmJh1nvzmPvUP4QuZBz\nTkT8oTzJ8hfgKvbOy/8HuJ28kNRpdzHt5thHryWdDVwJnBPlf3WUdSXPSvpRRNzeZYyq+cCt5COa\n7wA3kmsi6qYiiIiQdBE5crCBXIdwJfm/TIbjdZOE/Jh8sqdT7P+rJyJ+LulN4KfkOpH/AuvIhKyd\nNcA3Sj+NjYhFnRoUEc9JWgdcSC5ahUzq7i7762l+lPVjn72hz4YtoCS45O/dTexNmgAoic9s8qkX\nswOSIrq+CTGzzzlJ44BXgfkRcefnJXZLO2YBN0bE9NFqQx1JlwOzI+K8xsJmo8QjHmZWS9IpwEnk\n0xRHkYtnA7j/sxy7k4h4WNIUScd1WDsyWnaR02xmBywnHmbW5GpyEeMuYC1wVh8XLo5m7FoRsWy0\n29BORKwY7TaYNfFUi5mZmfWNH6c1MzOzvnHiYWZmZn3jxMPMzMz6xomHmZmZ9Y0TDzMzM+sbJx5m\nZmbWN048zMzMrG+ceJiZmVnfOPEwMzOzvnHiYWZmZn3zIbflgyhQ/Ap3AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f4594721d68>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_bs_decay()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\begin{pmatrix}\\mathbb{0}\\end{pmatrix}$"
      ],
      "text/plain": [
       "[[0]]"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "θ = Symbol('theta', real=True)\n",
    "slh.L.substitute({θ: 3*pi/2}).simplify_scalar()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "\\begin{align}\n",
       "  \\hat{H} &= \\hat{H}_0 + \\hat{H}_{\\text{int}} + \\hat{H}_{1} + \\hat{H}_{2} + \\hat{H}_{\\alpha}\\\\\n",
       "  \\hat{H}_0 &= - \\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} - \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})}\\\\\n",
       "  \\hat{H}_{\\text{int}} &= \\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{\\Pi}_{g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} + \\left(i \\kappa \\hat{a}^{(c_{1})\\dagger} \\otimes \\hat{a}^{(c_{2})} + \\text{c.c.}\\right) + \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{\\Pi}_{g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})}\\\\\n",
       "  \\hat{H}_{1} &= - \\frac{i \\Omega_{1} g_{1}}{2 \\Delta_{1}} \\left(\\hat{\\sigma}_{e,g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})} - \\text{c.c.}\\right)\\\\\n",
       "  \\hat{H}_{2} &= - \\frac{i \\Omega_{2} g_{2}}{2 \\Delta_{2}} \\left(\\hat{\\sigma}_{e,g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})} - \\text{c.c.}\\right)\\\\\n",
       "  \\hat{H}_{\\alpha} &= \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})} + \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})}\\\\\n",
       "\\end{align}"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "display_hamiltonian(slh.H)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "me = slh.symbolic_master_equation()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$- \\frac{1}{2} \\left(\\hat{\\rho}^{(q_{1} \\otimes c_{1} \\otimes q_{2} \\otimes c_{2})} \\otimes \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})\\dagger}\\right) \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})}\\right) + \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})\\dagger}\\right) \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})}\\right) \\otimes \\hat{\\rho}^{(q_{1} \\otimes c_{1} \\otimes q_{2} \\otimes c_{2})}\\right) - i \\left(-\\hat{\\rho}^{(q_{1} \\otimes c_{1} \\otimes q_{2} \\otimes c_{2})} \\left(\\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{\\Pi}_{g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} - \\frac{i \\Omega_{1} g_{1}}{2 \\Delta_{1}} \\hat{\\sigma}_{e,g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})} + \\frac{i \\Omega_{1} g_{1}}{2 \\Delta_{1}} \\hat{\\sigma}_{g,e}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})\\dagger} + \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})} - \\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} + i \\kappa \\hat{a}^{(c_{1})\\dagger} \\otimes \\hat{a}^{(c_{2})} - i \\kappa \\hat{a}^{(c_{1})} \\otimes \\hat{a}^{(c_{2})\\dagger} + \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{\\Pi}_{g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})} - \\frac{i \\Omega_{2} g_{2}}{2 \\Delta_{2}} \\hat{\\sigma}_{e,g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})} + \\frac{i \\Omega_{2} g_{2}}{2 \\Delta_{2}} \\hat{\\sigma}_{g,e}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})\\dagger} + \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})} - \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})}\\right) + \\left(\\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{\\Pi}_{g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} - \\frac{i \\Omega_{1} g_{1}}{2 \\Delta_{1}} \\hat{\\sigma}_{e,g}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})} + \\frac{i \\Omega_{1} g_{1}}{2 \\Delta_{1}} \\hat{\\sigma}_{g,e}^{(q_{1})} \\otimes \\hat{a}^{(c_{1})\\dagger} + \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{1})} - \\frac{g_{1}^{2}}{\\Delta_{1}} \\hat{a}^{(c_{1})\\dagger} \\hat{a}^{(c_{1})} + i \\kappa \\hat{a}^{(c_{1})\\dagger} \\otimes \\hat{a}^{(c_{2})} - i \\kappa \\hat{a}^{(c_{1})} \\otimes \\hat{a}^{(c_{2})\\dagger} + \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{\\Pi}_{g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})} - \\frac{i \\Omega_{2} g_{2}}{2 \\Delta_{2}} \\hat{\\sigma}_{e,g}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})} + \\frac{i \\Omega_{2} g_{2}}{2 \\Delta_{2}} \\hat{\\sigma}_{g,e}^{(q_{2})} \\otimes \\hat{a}^{(c_{2})\\dagger} + \\frac{\\sqrt{2} i \\lambda \\cos{\\left (\\theta \\right )} \\overline{\\Omega_{\\alpha}}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})\\dagger} - \\frac{\\sqrt{2} i \\Omega_{\\alpha} \\lambda \\cos{\\left (\\theta \\right )}}{2 \\sin{\\left (\\theta \\right )} - 2} \\hat{a}^{(c_{2})} - \\frac{g_{2}^{2}}{\\Delta_{2}} \\hat{a}^{(c_{2})\\dagger} \\hat{a}^{(c_{2})}\\right) \\hat{\\rho}^{(q_{1} \\otimes c_{1} \\otimes q_{2} \\otimes c_{2})}\\right) + \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})}\\right) \\otimes \\hat{\\rho}^{(q_{1} \\otimes c_{1} \\otimes q_{2} \\otimes c_{2})} \\otimes \\left(- \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{1})\\dagger} - \\frac{\\sqrt{2} \\sqrt{\\kappa} \\cos{\\left (\\theta \\right )}}{\\sin{\\left (\\theta \\right )} - 1} \\hat{a}^{(c_{2})\\dagger}\\right)$"
      ],
      "text/plain": [
       "-1/2 * (ρ̂^(q₁⊗c₁⊗q₂⊗c₂) ⊗ (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁)† - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)†) (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁) - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)) + (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁)† - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)†) (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁) - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)) ⊗ ρ̂^(q₁⊗c₁⊗q₂⊗c₂)) - ⅈ * (-ρ̂^(q₁⊗c₁⊗q₂⊗c₂) (g_1**2/Δ_1 * Π̂_g^(q₁) ⊗ â^(c₁)† â^(c₁) - I*Ω_1*g_1/(2*Δ_1) * σ̂_e,g^(q₁) ⊗ â^(c₁) + I*Ω_1*g_1/(2*Δ_1) * σ̂_g,e^(q₁) ⊗ â^(c₁)† + sqrt(2)*I*λ*cos(θ)*conjugate(Ω_α)/(2*(sin(θ) - 1)) * â^(c₁)† - sqrt(2)*I*Ω_α*λ*cos(θ)/(2*sin(θ) - 2) * â^(c₁) - g_1**2/Δ_1 * â^(c₁)† â^(c₁) + ⅈ⋅κ * â^(c₁)† ⊗ â^(c₂) - ⅈ⋅κ * â^(c₁) ⊗ â^(c₂)† + g_2**2/Δ_2 * Π̂_g^(q₂) ⊗ â^(c₂)† â^(c₂) - I*Ω_2*g_2/(2*Δ_2) * σ̂_e,g^(q₂) ⊗ â^(c₂) + I*Ω_2*g_2/(2*Δ_2) * σ̂_g,e^(q₂) ⊗ â^(c₂)† + sqrt(2)*I*λ*cos(θ)*conjugate(Ω_α)/(2*(sin(θ) - 1)) * â^(c₂)† - sqrt(2)*I*Ω_α*λ*cos(θ)/(2*sin(θ) - 2) * â^(c₂) - g_2**2/Δ_2 * â^(c₂)† â^(c₂)) + (g_1**2/Δ_1 * Π̂_g^(q₁) ⊗ â^(c₁)† â^(c₁) - I*Ω_1*g_1/(2*Δ_1) * σ̂_e,g^(q₁) ⊗ â^(c₁) + I*Ω_1*g_1/(2*Δ_1) * σ̂_g,e^(q₁) ⊗ â^(c₁)† + sqrt(2)*I*λ*cos(θ)*conjugate(Ω_α)/(2*(sin(θ) - 1)) * â^(c₁)† - sqrt(2)*I*Ω_α*λ*cos(θ)/(2*sin(θ) - 2) * â^(c₁) - g_1**2/Δ_1 * â^(c₁)† â^(c₁) + ⅈ⋅κ * â^(c₁)† ⊗ â^(c₂) - ⅈ⋅κ * â^(c₁) ⊗ â^(c₂)† + g_2**2/Δ_2 * Π̂_g^(q₂) ⊗ â^(c₂)† â^(c₂) - I*Ω_2*g_2/(2*Δ_2) * σ̂_e,g^(q₂) ⊗ â^(c₂) + I*Ω_2*g_2/(2*Δ_2) * σ̂_g,e^(q₂) ⊗ â^(c₂)† + sqrt(2)*I*λ*cos(θ)*conjugate(Ω_α)/(2*(sin(θ) - 1)) * â^(c₂)† - sqrt(2)*I*Ω_α*λ*cos(θ)/(2*sin(θ) - 2) * â^(c₂) - g_2**2/Δ_2 * â^(c₂)† â^(c₂)) ρ̂^(q₁⊗c₁⊗q₂⊗c₂)) + (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁) - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)) ⊗ ρ̂^(q₁⊗c₁⊗q₂⊗c₂) ⊗ (-sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₁)† - sqrt(2)*sqrt(κ)*cos(θ)/(sin(θ) - 1) * â^(c₂)†)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "me"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Numerical Values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "controls = sorted([sym for sym in slh.H.all_symbols() if sym.name.startswith('Omega')], key=str)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# from algebra import generate_num_vals_code\n",
    "# generate_num_vals_code(slh.all_symbols())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Delta =  5000.0  # MHz\n",
    "g     =    50.0  # MHz\n",
    "kappa =     0.5  # MHz\n",
    "E0    =    70.0  # MHz\n",
    "t0    = -4    # microsec\n",
    "T     =  4    # microsec\n",
    "nt    = 2001\n",
    "theta = 1.5 * np.pi\n",
    "num_vals = {\n",
    "    Symbol('Delta_1', real=True):   Delta,\n",
    "    Symbol('Delta_2', real=True):   Delta,\n",
    "    Symbol('g_2', positive=True):   g,\n",
    "    Symbol('g_1', positive=True):   g,\n",
    "    Symbol('kappa', positive=True): kappa,\n",
    "    Symbol('theta', real=True):     theta,\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Magnitude of coefficients"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from qnet.algebra.pattern_matching import pattern\n",
    "from qnet.algebra.operator_algebra import ScalarTimesOperator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'-I*Omega_1*g_1/(2*Delta_1)': -0.005*I*Omega_1,\n",
       " '-I*Omega_2*g_2/(2*Delta_2)': -0.005*I*Omega_2,\n",
       " '-I*kappa': -0.5*I,\n",
       " '-g_1**2/Delta_1': -0.500000000000000,\n",
       " '-g_2**2/Delta_2': -0.500000000000000,\n",
       " '-sqrt(2)*I*Omega_alpha*lambda*cos(theta)/(2*sin(theta) - 2)': -4.59242549680257e-17*sqrt(2)*I*Omega_alpha*lambda,\n",
       " 'I*Omega_1*g_1/(2*Delta_1)': 0.005*I*Omega_1,\n",
       " 'I*Omega_2*g_2/(2*Delta_2)': 0.005*I*Omega_2,\n",
       " 'I*kappa': 0.5*I,\n",
       " 'g_1**2/Delta_1': 0.500000000000000,\n",
       " 'g_2**2/Delta_2': 0.500000000000000,\n",
       " 'sqrt(2)*I*lambda*cos(theta)*conjugate(Omega_alpha)/(2*(sin(theta) - 1))': 4.59242549680257e-17*sqrt(2)*I*lambda*conjugate(Omega_alpha)}"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "{str(t.coeff): t.coeff.subs(num_vals) for t in pattern(head=ScalarTimesOperator).findall(slh.H)}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Time Grid and Guess Pulse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "tgrid = QDYN.pulse.pulse_tgrid(t0=t0, T=T, nt=nt) # microsec"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "p = QDYN.pulse.Pulse(\n",
    "    tgrid, amplitude=(E0 * blackman(tgrid, t0, T)),\n",
    "    time_unit='microsec', ampl_unit='MHz')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pa = QDYN.pulse.Pulse(\n",
    "    tgrid, amplitude=(0.0 * blackman(tgrid, t0, T)),\n",
    "    time_unit='microsec', ampl_unit='MHz')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "controls = OrderedDict([\n",
    "    (Symbol('Omega_1'), p.copy()),\n",
    "    (Symbol('Omega_2'), p.copy()),\n",
    "    (Symbol('Omega_a'), pa.copy()),\n",
    "])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## OCT for forward transfer (hs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "Scalar coefficient '4.59242549680257e-17*sqrt(2)*I*lambda*conjugate(Omega_alpha)' is not numerical",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m/home/goerz/Documents/Programming/QNET/qnet/convert/to_qutip.py\u001b[0m in \u001b[0;36mconvert_to_qutip\u001b[0;34m(expr, full_space, mapping)\u001b[0m\n\u001b[1;32m    122\u001b[0m         \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 123\u001b[0;31m             \u001b[0mcoeff\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcomplex\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcoeff\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    124\u001b[0m         \u001b[0;32mexcept\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/.conda/envs/default/lib/python3.5/site-packages/sympy/core/expr.py\u001b[0m in \u001b[0;36m__complex__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    230\u001b[0m         \u001b[0mre\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mim\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mresult\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mas_real_imag\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 231\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0mcomplex\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfloat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mre\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfloat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mim\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    232\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/.conda/envs/default/lib/python3.5/site-packages/sympy/core/expr.py\u001b[0m in \u001b[0;36m__float__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    225\u001b[0m             \u001b[0;32mraise\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"can't convert complex to float\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 226\u001b[0;31m         \u001b[0;32mraise\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"can't convert expression to float\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    227\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mTypeError\u001b[0m: can't convert expression to float",
      "\nDuring handling of the above exception, another exception occurred:\n",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-25-c2272798f918>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      2\u001b[0m     \u001b[0mslh\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnum_vals\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcontrols\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0menergy_unit\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'MHz'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m     \u001b[0mmcwf\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mFalse\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnon_herm\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mTrue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0moct_target\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'excitation_transfer_fw'\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 4\u001b[0;31m     lambda_a=1e-4)\n\u001b[0m",
      "\u001b[0;32m/home/goerz/quantum_transmission2/src/qdyn_model_v1.py\u001b[0m in \u001b[0;36mmake_qdyn_oct_model\u001b[0;34m(network_slh, num_vals, controls, oct_target, energy_unit, mcwf, non_herm, use_dissipation_superop, lambda_a, seed, allow_negative, **kwargs)\u001b[0m\n\u001b[1;32m    393\u001b[0m         \u001b[0mnetwork_slh\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnum_vals\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcontrols\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0menergy_unit\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0menergy_unit\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmcwf\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mmcwf\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    394\u001b[0m         \u001b[0mnon_herm\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mnon_herm\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0muse_dissipation_superop\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0muse_dissipation_superop\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 395\u001b[0;31m         states=states)\n\u001b[0m\u001b[1;32m    396\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    397\u001b[0m     \u001b[0;32mfor\u001b[0m \u001b[0mi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpulse\u001b[0m \u001b[0;32min\u001b[0m \u001b[0menumerate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpulses\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/quantum_transmission2/src/qdyn_model_v1.py\u001b[0m in \u001b[0;36mmake_qdyn_model\u001b[0;34m(network_slh, num_vals, controls, energy_unit, mcwf, non_herm, use_dissipation_superop, states, nodiss)\u001b[0m\n\u001b[1;32m    234\u001b[0m     \u001b[0mn_nodes\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0ml\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ml\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mhs_labels\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0ml\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstartswith\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'q'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    235\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 236\u001b[0;31m     H0 = (convert_to_qutip(ham_parts['H0'], full_space=hs) +\n\u001b[0m\u001b[1;32m    237\u001b[0m           convert_to_qutip(ham_parts['Hint'], full_space=hs))\n\u001b[1;32m    238\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/Documents/Programming/QNET/qnet/convert/to_qutip.py\u001b[0m in \u001b[0;36mconvert_to_qutip\u001b[0;34m(expr, full_space, mapping)\u001b[0m\n\u001b[1;32m    118\u001b[0m         \u001b[0;32mreturn\u001b[0m \u001b[0m_convert_local_operator_to_qutip\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfull_space\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmapping\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    119\u001b[0m     \u001b[0;32melif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mOperatorOperation\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 120\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0m_convert_operator_operation_to_qutip\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfull_space\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmapping\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    121\u001b[0m     \u001b[0;32melif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mScalarTimesOperator\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    122\u001b[0m         \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/Documents/Programming/QNET/qnet/convert/to_qutip.py\u001b[0m in \u001b[0;36m_convert_operator_operation_to_qutip\u001b[0;34m(expr, full_space, mapping)\u001b[0m\n\u001b[1;32m    264\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mOperatorPlus\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    265\u001b[0m         return sum((convert_to_qutip(op, full_space, mapping=mapping)\n\u001b[0;32m--> 266\u001b[0;31m                     for op in expr.operands), 0)\n\u001b[0m\u001b[1;32m    267\u001b[0m     \u001b[0;32melif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mOperatorTimes\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    268\u001b[0m         \u001b[0;31m# if any factor acts non-locally, we need to expand distributively.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/Documents/Programming/QNET/qnet/convert/to_qutip.py\u001b[0m in \u001b[0;36m<genexpr>\u001b[0;34m(.0)\u001b[0m\n\u001b[1;32m    264\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mOperatorPlus\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    265\u001b[0m         return sum((convert_to_qutip(op, full_space, mapping=mapping)\n\u001b[0;32m--> 266\u001b[0;31m                     for op in expr.operands), 0)\n\u001b[0m\u001b[1;32m    267\u001b[0m     \u001b[0;32melif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexpr\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mOperatorTimes\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    268\u001b[0m         \u001b[0;31m# if any factor acts non-locally, we need to expand distributively.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/goerz/Documents/Programming/QNET/qnet/convert/to_qutip.py\u001b[0m in \u001b[0;36mconvert_to_qutip\u001b[0;34m(expr, full_space, mapping)\u001b[0m\n\u001b[1;32m    124\u001b[0m         \u001b[0;32mexcept\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    125\u001b[0m             raise TypeError(\"Scalar coefficient '%s' is not numerical\" %\n\u001b[0;32m--> 126\u001b[0;31m                             expr.coeff)\n\u001b[0m\u001b[1;32m    127\u001b[0m         return coeff * convert_to_qutip(expr.term, full_space=full_space,\n\u001b[1;32m    128\u001b[0m                                         mapping=mapping)\n",
      "\u001b[0;31mTypeError\u001b[0m: Scalar coefficient '4.59242549680257e-17*sqrt(2)*I*lambda*conjugate(Omega_alpha)' is not numerical"
     ]
    }
   ],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals, controls, energy_unit='MHz',\n",
    "    mcwf=False, non_herm=True, oct_target='excitation_transfer_fw',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data_root = \"data/DR_BS_FB/\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(data_root+\"fw_oct_hs001/\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for $\\sqrt{SWAP}$ (hs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sqrt_SWAP = QDYN.gate2q.Gate2Q([\n",
    "    [1,             0,             0,   0],\n",
    "    [0, -1/np.sqrt(2),  1/np.sqrt(2),   0],\n",
    "    [0, -1/np.sqrt(2), -1/np.sqrt(2),   0],\n",
    "    [0,             0,             0, -1j]]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals, controls, energy_unit='MHz',\n",
    "    mcwf=False, non_herm=True, oct_target='gate',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/sqrtSWAP_oct_hs001/\")\n",
    "sqrt_SWAP.write(\"./DR_BS_FB/sqrtSWAP_oct_hs001/target_gate.dat\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for SWAP (hs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals, controls, energy_unit='MHz',\n",
    "    mcwf=False, non_herm=True, oct_target='gate',\n",
    "    lambda_a=1e-2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/SWAP_oct_hs001/\")\n",
    "QDYN.gate2q.SWAP.write(\"./DR_BS_FB/SWAP_oct_hs001/target_gate.dat\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for PE (hs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals, controls, energy_unit='MHz',\n",
    "    mcwf=False, non_herm=True, oct_target='gate',\n",
    "    lambda_a=1e-2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/PE_oct_hs001/\")\n",
    "QDYN.gate2q.CPHASE.write(\"./DR_BS_FB/PE_oct_hs001/target_gate.dat\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for sqrtSWAP (LI, hs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals, controls, energy_unit='MHz',\n",
    "    mcwf=False, non_herm=True, oct_target='gate',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/sqrtSWAP_oct_hsLI001/\")\n",
    "QDYN.gate2q.sqrt_SWAP.write(\"./DR_BS_FB/sqrtSWAP_oct_hsLI001/target_gate.dat\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for sqrtSWAP (hs, dissipative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "slh = single_sided_network.network_slh(n_cavity=2, n_nodes=2, topology='bs_fb')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sqrt_SWAP = QDYN.gate2q.Gate2Q([\n",
    "    [1,             0,             0,   0],\n",
    "    [0, -1/np.sqrt(2),  1/np.sqrt(2),   0],\n",
    "    [0, -1/np.sqrt(2), -1/np.sqrt(2),   0],\n",
    "    [0,             0,             0, -1j]]\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "num_vals_diss = num_vals.copy()\n",
    "θ = Symbol('theta', real=True)\n",
    "num_vals_diss[θ] = 0.25 * np.pi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals_diss, controls, energy_unit='MHz',\n",
    "    mcwf=True, non_herm=True, oct_target='gate',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/sqrtSWAP_diss_oct_hs001/\")\n",
    "sqrt_SWAP.write(\"./DR_BS_FB/sqrtSWAP_diss_oct_hs001/target_gate.dat\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for forward transfer (hs, dissipative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals_diss, controls, energy_unit='MHz',\n",
    "    mcwf=True, non_herm=True, oct_target='excitation_transfer_fw',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "L = model.lindblad_ops()[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/fw_diss_oct_hs001/\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  OCT for dicke state initialization (hs, dissipative)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model = qdyn_model.make_qdyn_oct_model(\n",
    "    slh, num_vals_diss, controls, energy_unit='MHz',\n",
    "    mcwf=True, non_herm=True, oct_target='dicke',\n",
    "    lambda_a=1e-4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model.write_to_runfolder(\"./DR_BS_FB/dicke_diss_oct_hs001/\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example optimization commands\n",
    "\n",
    "*   hilbert space optimization\n",
    "\n",
    "    ```OMP_NUM_THREADS=1 qdyn_optimize --internal-units=GHz_units.txt --debug --J_T=J_T_re .```\n",
    "\n",
    "*   trajectory optimization (model with `mcwf=True`)\n",
    "\n",
    "    ```OMP_NUM_THREADS=1 qdyn_optimize --internal-units=GHz_units.txt --n-trajs=20 --J_T=J_T_re .```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [default]",
   "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.5.3"
  },
  "toc": {
   "nav_menu": {
    "height": "123px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": 4,
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}