{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Incoherent Type 4 FFL with AND gate\n", "\n", "Implement the incoherent type 4 FFL with an AND gate, as described by Mangan & Alon (PNAS, 2003).\n", "\n", "The parameters are\n", "- B_y: basal rate of Y's transcription\n", "- Beta_y: maximum rate of Y's transcription\n", "- B_z: basal rate of Z's transcription\n", "- Beta_z: maximum rate of Z's transcription\n", "- K_xy: the activation threshold for X's activation of Y\n", "- K_xz: the activation threshold for X's activation of Z\n", "- K_yz: the activation threshold for Y's activation of Z\n", "- H: the Hill coefficient (exponent) for all regulation\n", "- alpha_y: the degradation rate constant for Y\n", "- alpha_z: the degradation rate constant for Z\n", "- SX: amount of active X in system (X*=SX)\n", "- SY: whether or not Y is activated (Y*=Y if SY=1, Y*=0 if SY=0)\n", "\n", "The model is\n", "\\begin{eqnarray*}\n", "\\frac{dY}{dt} &=& B_y + \\beta_Y \\frac{1}{1+(X^*/K_{xy})^H} - \\alpha_y Y\\\\\n", "\\frac{dZ}{dt} &=& B_z + \\beta_Z \\left(\\frac{(X^*/K_{xz})^H}{1+(X^*/K_{xz})^H}\\right) \\left(\\frac{(Y^*/K_{yz})^H}{1+(Y^*/K_{yz})^H}\\right) - \\alpha_z Z\\\\\n", "\\end{eqnarray*}\n" ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "def incoherent_type_4_and_model( t, statevars, params ):\n", " Y = statevars[0];\n", " Z = statevars[1];\n", " B_y = params[0]\n", " Beta_y = params[1]\n", " B_z = params[2]\n", " Beta_z = params[3]\n", " K_xy = params[4]\n", " K_xz = params[5]\n", " K_yz = params[6]\n", " H = params[7]\n", " alpha_y = params[8]\n", " alpha_z = params[9]\n", " SX = params[10]\n", " SY = params[11]\n", " \n", " if SX == 0:\n", " Xstar = 0\n", " else:\n", " Xstar = 1\n", " if SY == 0:\n", " Ystar = 0\n", " else:\n", " Ystar = Y\n", " \n", " dydt = np.zeros( (2,) )\n", " yreg = 1/(1 + (Xstar/K_xy)**H)\n", " zregx = ((Xstar/K_xz)**H)/(1 + (Xstar/K_xz)**H)\n", " zregy = ((Ystar/K_yz)**H)/(1 + (Ystar/K_yz)**H)\n", " zreg = zregx * zregy\n", " dydt[0] = B_y + Beta_y*yreg - alpha_y*Y;\n", " dydt[1] = B_z + Beta_z*zreg - alpha_z*Z;\n", "\n", " return dydt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Simulate the model\n", "\n", "Simulate the model. In the paper, they specify that for most simulations,\n", "\n", "- B_y = 0\n", "- Beta_y = 1\n", "- B_z = 0\n", "- Beta_z = 1\n", "- H = 2\n", "- alpha_y = 1\n", "- alpha_z = 1\n", "\n", "The activation thresholds are given in the figure captions, as are SX and SY." ] }, { "cell_type": "code", "execution_count": 38, "metadata": {}, "outputs": [], "source": [ "import scipy.integrate\n", "\n", "B_y = 0\n", "Beta_y = 1\n", "B_z = 0\n", "Beta_z = 1\n", "K_xy = 0.1\n", "K_xz = 0.1\n", "K_yz = 1 # 1, 0.3, 0.1\n", "H = 2\n", "alpha_y = 1\n", "alpha_z = 1\n", "\n", "time_x_on = 0.5\n", "time_x_off = 2.5\n", "params_off = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,0,1);\n", "params_on = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,1,1);\n", "\n", "mf_off = lambda t,y : incoherent_type_4_and_model(t,y, params_off)\n", "mf_on = lambda t,y : incoherent_type_4_and_model(t,y, params_on)\n", "\n", "sol1 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(0,time_x_on), y0=(0,0), \n", " t_eval=np.arange(0,time_x_on,.01) )\n", "sol2 = scipy.integrate.solve_ivp( fun=mf_on, t_span=(time_x_on,time_x_off), y0=sol1.y[:,-1], \n", " t_eval=np.arange(time_x_on,time_x_off,.01) )\n", "sol3 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(time_x_off,6), y0=sol2.y[:,-1], \n", " t_eval=np.arange(time_x_off,6,.01) )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plot the simulation\n", "Plot X over time in an upper subplot and Z over time in a lower subplot." ] }, { "cell_type": "code", "execution_count": 39, "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAAArlElEQVR4nO3deZgV5Zn38e/dOzR009DN2kAjIIiooC3K4hIdHTQGMokT0SzGMcEkmmjyzmQ0M4mJc+WdmHcmZtFkJErijgY3MiFuQSMqII2ibAItgnQLdLPTNPR6v3+cgjQtq56i6pz+fa7rXOfUU1Xn3HWh53eq6unnMXdHREQkbjKiLkBERORgFFAiIhJLCigREYklBZSIiMSSAkpERGIpK+oCkqW4uNjLysqiLkNERI7RokWLNrt7Sfv2tAmosrIyKioqoi5DRESOkZmtO1h7aJf4zGy6mdWY2dJDrDcz+6WZVZrZ22Z2ept1V5vZ6uBxdVg1iohIfIV5D+r3wMTDrL8EGBo8pgK/ATCz7sCtwFnAGOBWMysKsU4REYmh0C7xufvLZlZ2mE0mA/d7YiiL+WbWzcz6AOcDz7v7VgAze55E0D0SVq1ReHjB+zy9uDrqMuQjmjyqH1edNSDqMkTSWpS9+PoB69ssVwVth2r/EDObamYVZlZRW1sbWqFheHpxNcs37Iy6DPkIlm/YqR8XIsdBSneScPdpwDSA8vLylBtUcESfAh69bmzUZcgxuuLueVGXINIhRHkGVQ30b7NcGrQdql1ERDqQKANqFvCloDff2cAOd98APAtcbGZFQeeIi4M2ERHpQEK7xGdmj5Do8FBsZlUkeuZlA7j7/wCzgUuBSqAeuCZYt9XM/gNYGLzVbfs6TIiISMcRZi++K4+w3oHrD7FuOjA9jLpERCQ1aCw+ERGJJQWUiIjEkgJKRERiSQElIiKxpIASEZFYUkCJiEgsKaBERCSWFFAiIhJLCigREYklBZSIiMSSAkpERGJJASUiIrGkgBIRkVhSQImISCwpoEREJJZCDSgzm2hmK82s0sxuPsj6O8xscfBYZWbb26xrabNuVph1iohI/IQ5o24mcBdwEVAFLDSzWe6+fN827v7tNtt/Exjd5i32uPuosOoTEZF4C/MMagxQ6e5r3L0RmAFMPsz2VwKPhFiPiIikkDADqh+wvs1yVdD2IWY2EBgEzGnTnGdmFWY238w+fYj9pgbbVNTW1iapbBERiYO4dJKYAsx095Y2bQPdvRy4Cvi5mQ1uv5O7T3P3cncvLykpOV61iojIcRBmQFUD/dsslwZtBzOFdpf33L06eF4DvMSB96dERCTNhRlQC4GhZjbIzHJIhNCHeuOZ2XCgCJjXpq3IzHKD18XAeGB5+31FRCR9hdaLz92bzewG4FkgE5ju7svM7Dagwt33hdUUYIa7e5vdTwLuNrNWEiH6k7a9/0REJP2FFlAA7j4bmN2u7Qftln94kP1eA04JszYREYm3uHSSEBEROYACSkREYkkBJSIisaSAEhGRWFJAiYhILCmgREQklhRQIiISSwooERGJJQWUiIjEkgJKRERiSQElIiKxpIASEZFYUkCJiEgsKaBERCSWQg0oM5toZivNrNLMbj7I+i+bWa2ZLQ4eX2mz7mozWx08rg6zThERiZ/Q5oMys0zgLuAioApYaGazDjLx4KPufkO7fbsDtwLlgAOLgn23hVWviIjEyyHPoMxswGHWnXMU7z0GqHT3Ne7eCMwAJh9lXX8PPO/uW4NQeh6YeJT7iohIGjjcJb6XzOy7wZkQAGbWy8weBO44ivfuB6xvs1wVtLX3WTN728xmmln/Y9nXzKaaWYWZVdTW1h5FSSIikioOF1BnAIOBxWZ2gZndCLwOzCNxdpQMfwTK3P1UEmdJ9x3Lzu4+zd3L3b28pKQkSSWJiEgcHPIeVHBp7bogmF4APgDOdveqo3zvaqB/m+XSoK3tZ2xps3gP8NM2+57fbt+XjvJzRUQkDRzuHlQ3M7sbuIbE/Z+ZwJ/N7IKjfO+FwFAzG2RmOcAUYFa7z+jTZnESsCJ4/SxwsZkVmVkRcHHQJiIiHcThevG9AfwauN7dm4HnzGwU8GszW+fuVx7ujd292cxuIBEsmcB0d19mZrcBFe4+C/iWmU0CmoGtwJeDfbea2X+QCDmA29x960c+ShERSTmHC6hz21/Oc/fFwDgz++rRvLm7zwZmt2v7QZvXtwC3HGLf6cD0o/kcERFJP4e8xHe4e03u/ttwyhEREUnQUEciIhJLCigREYklBZSIiMSSAkpERGJJASUiIrGkgBIRkVhSQImISCwpoEREJJYUUCIiEksKKBERiSUFlIiIxJICSkREYkkBJSIisaSAEhGRWAo1oMxsopmtNLNKM7v5IOu/Y2bLzextM/uLmQ1ss67FzBYHj1nt9xURkfR2uAkLPxYzywTuAi4CqoCFZjbL3Ze32exNoNzd683s68BPgSuCdXvcfVRY9YmISLyFeQY1Bqh09zXu3gjMACa33cDdX3T3+mBxPlAaYj0iIpJCwgyofsD6NstVQduhXAv8uc1ynplVmNl8M/v0wXYws6nBNhW1tbUfu2AREYmP0C7xHQsz+wJQDpzXpnmgu1eb2QnAHDNb4u7vtt3P3acB0wDKy8v9uBUsIiKhC/MMqhro32a5NGg7gJn9HfBvwCR3b9jX7u7VwfMa4CVgdIi1iohIzIQZUAuBoWY2yMxygCnAAb3xzGw0cDeJcKpp015kZrnB62JgPNC2c4WIiKS50C7xuXuzmd0APAtkAtPdfZmZ3QZUuPss4P8BXYA/mBnA++4+CTgJuNvMWkmE6E/a9f4TEZE0F+o9KHefDcxu1/aDNq//7hD7vQacEmZtIiISbxpJQkREYkkBJSIisaSAEhGRWFJAiYhILCmgREQklhRQIiISSwooERGJJQWUiIjEkgJKRERiSQElIiKxpIASEZFYUkCJiEgsKaBERCSWFFAiIhJLoQaUmU00s5VmVmlmNx9kfa6ZPRqsX2BmZW3W3RK0rzSzvw+zThERiZ/QAsrMMoG7gEuAEcCVZjai3WbXAtvcfQhwB3B7sO8IEjPwngxMBH4dvJ+IiHQQYU5YOAaodPc1AGY2A5jMgVO3TwZ+GLyeCdxpial1JwMz3L0BeM/MKoP3mxdWsT/64zKWf7AzrLf/kOUbdjKiT8Fx+zxJruUbdnLF3aH95yiSMkb0LeDWT50cynuHGVD9gPVtlquAsw61TTBF/A6gR9A+v92+/dp/gJlNBaYCDBgwIGmFHw8j+hQwedSHDklSgP7dRI6PUKd8D5u7TwOmAZSXl/vHea+wfgFI+rnqrAFcdVZq/SASSUVhdpKoBvq3WS4N2g66jZllAYXAlqPcV0RE0liYAbUQGGpmg8wsh0Snh1nttpkFXB28vhyY4+4etE8JevkNAoYCr4dYq4iIxExol/iCe0o3AM8CmcB0d19mZrcBFe4+C7gXeCDoBLGVRIgRbPcYiQ4VzcD17t4SVq0iIhI/ljhhSX1mVgus+5hvUwxsTkI5caXjS13pfGyg40t1H/f4Brp7SfvGtAmoZDCzCncvj7qOsOj4Ulc6Hxvo+FJdWMenoY5ERCSWFFAiIhJLCqgDTYu6gJDp+FJXOh8b6PhSXSjHp3tQIiISSzqDEhGRWFJAiYhILCmgAkeauyqVmdl0M6sxs6VR15JsZtbfzF40s+VmtszMboy6pmQyszwze93M3gqO70dR15RsZpZpZm+a2f9GXUuymdlaM1tiZovNrCLqepLNzLqZ2Uwze8fMVpjZ2KS+v+5B7Z+7ahVwEYmR0xcCV7r78sPumCLM7FygDrjf3UdGXU8ymVkfoI+7v2FmXYFFwKfT6N/OgHx3rzOzbOAV4EZ3n3+EXVOGmX0HKAcK3P2yqOtJJjNbC5S7e1r+ka6Z3QfMdfd7giHtOrv79mS9v86gEvbPXeXujcC+uavSgru/TGIoqbTj7hvc/Y3g9S5gBQeZmiVVeUJdsJgdPNLmV6WZlQKfBO6JuhY5NmZWCJxLYsg63L0xmeEECqh9DjZ3Vdp8yXUUZlYGjAYWRFxKUgWXwBYDNcDz7p5Ox/dz4LtAa8R1hMWB58xsUTB/XToZBNQCvwsu0d5jZvnJ/AAFlKQFM+sCPA7c5O7Hb2rk48DdW9x9FIlpZ8aYWVpcpjWzy4Aad18UdS0hmuDupwOXANcHl9vTRRZwOvAbdx8N7AaSev9eAZWg+adSWHBv5nHgIXd/Iup6whJcPnkRmBhxKckyHpgU3KeZAVxgZg9GW1JyuXt18FwDPEnidkK6qAKq2pzRzyQRWEmjgEo4mrmrJIaCTgT3Aivc/WdR15NsZlZiZt2C151IdOR5J9KiksTdb3H3UncvI/H/3Bx3/0LEZSWNmeUHHXcILn1dDKRNT1p33wisN7NhQdOFJKZISpqUnvI9WQ41d1XEZSWNmT0CnA8Um1kVcKu73xttVUkzHvgisCS4TwPwPXefHV1JSdUHuC/oaZoBPObuadcdO031Ap5M/IYiC3jY3Z+JtqSk+ybwUPDDfg1wTTLfXN3MRUQklnSJT0REYkkBJSIisaSAEhGRWFJAiYhILCmgREQklhRQIiISSwooERGJJQWUiIjEkgJKRERiSQElIiKxpIASEZFYSpvBYouLi72srCzqMkRE5BgtWrRos7uXtG9Pm4AqKyujoqIi6jJEROQYmdm6g7XrEp+IiMSSAkoOS9OxiEhU0uYSn3w8lTW7+MuKGhau3ca7tXVUb99DY3MrAF1zs+hdmEdZcT6j+nfjjIFFnDGwiOxM/b4RkfAooDqw5pZW/rRkA/fMfY8l1TsAOKE4nxF9C7hoRC/ysjMxYMeeJjbs2MPqmjqeX74JgIK8LC48qRefPb2UcYN7kJFhER6JiKQjBVQHNe/dLXz/6aVU1tQxtGcXfnDZCD55ah96FeQddr8d9U3MW7OF55dv4oUVm3jyzWrKenTmqrMGcOWYAXTNyz5ORyAi6S5tpnwvLy939eI7sobmFn78pxXcP28d/bt34t8uPYmLR/T+SGdAe5taeGbpRh5e8D6vr91KYadsrhlfxjXjBlHYWUElIkfHzBa5e/mH2hVQHUfNzr189f4K3qrawbUTBvHPFw+jU05mUt777art/GpOJc8v30Rhp2xuvHAoXxw7UPepROSIFFAd3Lotu/nCvQvYUtfIzz43iokje4fyOcs/2Ml//nkFc1dv5oTifL536UlceFJPzHSPSkQO7lABpZ+3HcD6rfV87u551O1t5uGvnh1aOAGM6FvA/f80hulfLgeDr9xfwXUPLGLjjr2hfaaIpCcFVJrbUtfA1dNfZ09jC49MPZtR/buF/plmxgXDe/HsTedy8yXD+euqWi762V95aME6WlvT44xdRMKngEpjTS2tfP3BN6jevofpXz6T4b0LjuvnZ2dm8LXzBvPsTecysl8h//bkUqb8dj7rt9Yf1zpEJDUpoNLYf85+h9fXbuWnl59KeVn3yOooK87n4a+exU8/eyorPtjJJb+Yy+OLqjRKhYgclgIqTT2/fBPTX32PL48rY/KoflGXg5nxuTP7M/vGczipT1f+zx/e4oZH3mRHfVPUpYlITIUaUGY20cxWmlmlmd18kPW5ZvZosH6BmZUF7WPMbHHweMvM/iHMOtPN1t2N3PLE25zUp4DvXXpS1OUcoH/3zsyYOpZ/+fthPLt0IxN/8TIL1myJuiwRiaHQAsrMMoG7gEuAEcCVZjai3WbXAtvcfQhwB3B70L4UKHf3UcBE4G4z06gXR+kHTy9lx54mfva508jJit9JcmaGcf0nhvD418eRl53Jlb+dz29eelcdKETkAGF+e40BKt19jbs3AjOAye22mQzcF7yeCVxoZubu9e7eHLTnAfrmOkovr6rlf9/ewDcvGMpJfY5vp4hjdVr/bsy6YTyXjOzD7c+8w1fvr9AlPxHZL8yA6gesb7NcFbQddJsgkHYAPQDM7CwzWwYsAb7WJrD2M7OpZlZhZhW1tbUhHEJqaWxu5Yd/XEZZj85cd94JUZdzVLrmZXPnVaO59VMjeHl1LZ/81VzertoedVkiEgPxu/4TcPcF7n4ycCZwi5l9aBRTd5/m7uXuXl5S8qHZgjuc+15by5ra3dz6qZPJzUrOEEbHg5lxzfhBPHbdWFpbnct/M48H5q9TLz+RDi7MgKoG+rdZLg3aDrpNcI+pEDjgjrm7rwDqgJGhVZoGdu5t4q6XKjn3xBI+Mbxn1OV8JKMHFPGnb53DuCE9+P5TS/nuzLfZ29QSdVkiEpEwA2ohMNTMBplZDjAFmNVum1nA1cHry4E57u7BPlkAZjYQGA6sDbHWlHfP3PfYXt/Ev1w8LOpSPpai/BymX30mN144lD8squKKafM1TJJIBxVaQAX3jG4AngVWAI+5+zIzu83MJgWb3Qv0MLNK4DvAvq7oE4C3zGwx8CTwDXffHFatqW7r7kbunbuGS0/pzSmlhVGX87FlZBjfvuhE7v7iGVRu2sVlv3qFirVboy5LRI4zjWaeBn72/Cp++ZfVPP/tcxnaq2vU5STVqk27mHp/BdXb9/DDSSfz+bMGRl2SiCSZRjNPU/WNzdw/by0XjeiVduEEcGKvrjx9/QTGDynm355cyi1PLKGhWfelRDoCBVSKe3TherbXN/G18wZHXUpoCjtnc+/VZ/KN8wfzyOvvc9VvF1CzU/elRNKdAiqFtbQ698x9jzPLijhjYFHU5YQqM8P47sTh3HnVaJZ/sJNP3fkKb76/LeqyRCRECqgU9tdVNVRv38M/jR8UdSnHzWWn9uWJb4wjJyuDK+6ez2MV64+8k4ikJAVUCnto/vuUdM3l70b0irqU4+qkPgXMun4CZw4q4rsz3+aHs5bR1NIadVkikmQKqBRVta2eOStruKK8P9mZHe+fsSg/h/uuGcNXJgzi96+t5Qv3LGBLXUPUZYlIEnW8b7Y08ejCxKWtKWP6H2HL9JWVmcG/XzaCO644jcXrtzPpzldZWr0j6rJEJEkUUCmoqaWVRxeu5xPDelJa1DnqciL3D6NLmfm1cbg7l//Pazy9uP2IWiKSihRQKeiV1Zup2dXAlDM77tlTe6eUFjLrmxM4tV83bpyxmP87ewXNui8lktIUUCno6cXVFHbK5vxhqTkobFiKu+Ty4FfO4ktjBzLt5TVc8/uFbK9vjLosEfmIFFAppr6xmeeWb+LSU/rEcrbcqOVkZXDb5JH85DOnMH/NFibd+SrvbNwZdVki8hHoGy7FvLCihvrGFiaP6ht1KbE2ZcwAZkwdy96mFj7z69f485INUZckIsdIAZViZi2upndBHmPKukddSuydMbCIP35zAif26srXH3qD/35uJa2t6TE4skhHoIBKIdt2N/LSylomjepLRoZFXU5K6FWQx6PXnc3nykv51ZxKvnJ/he5LiaQIBVQKeW75RppbnUmn6fLescjNyuT2z57KbZNPZu7qWi771Su8XbU96rJE5AgUUCnk+eWb6NetEyf3LYi6lJRjZnxpbBmPXTeW1lbn8t/M44H560iX+dBE0pECKkXUNzYzd/VmLhrRCzNd3vuoRg8o4k/fOoexg3vw/aeW8u1HF7O7oTnqskTkIBRQKWLu6s00NLdyUQcbGDYMRfk5/O7LZ/LPF5/IrLc+YPJdr1JZsyvqskSknVADyswmmtlKM6s0s5sPsj7XzB4N1i8ws7Kg/SIzW2RmS4LnC8KsMxW8sHwTXfOyGDNIvfeSISPDuOGCoTxw7Vlsr29k0p2vaogkkZgJLaDMLBO4C7gEGAFcaWYj2m12LbDN3YcAdwC3B+2bgU+5+ynA1cADYdWZClpanTnv1PCJYT075MjlYRo/pJg/fescTu5bwI0zFvO9J5ewt0lTyovEQZjfdmOASndf4+6NwAxgcrttJgP3Ba9nAheambn7m+7+QdC+DOhkZrkh1hprb76/jS27G3V5LyS9CvJ4+Ktnc915J/DwgveZdOcrGn1CJAbCDKh+QNvpTquCtoNu4+7NwA6gR7ttPgu84e4fmuzHzKaaWYWZVdTW1iat8Lh5fsUmsjON84aVRF1K2srOzOCWS07i/n8aw9bdTUy681Xun7dWvfxEIhTr60VmdjKJy37XHWy9u09z93J3Ly8pSd8v77mrNnPGwCIK8rKjLiXtnXtiCc/cdA7jBvfgB08v46v3L2Lrbv1hr0gUwgyoaqDtfBClQdtBtzGzLKAQ2BIslwJPAl9y93dDrDPWNtc1sHzDTs4Zmr4BHDfFXXKZfvWZfP+yEby8qpZLfvEyr727OeqyRDqcMANqITDUzAaZWQ4wBZjVbptZJDpBAFwOzHF3N7NuwJ+Am9391RBrjL1XKxNfjBOGFEdcSceSkWFcO2EQT3xjHPm5WXz+ngX85+wVNDSrA4XI8RJaQAX3lG4AngVWAI+5+zIzu83MJgWb3Qv0MLNK4DvAvq7oNwBDgB+Y2eLg0SEnP3pl9WYKO2Uzsl9h1KV0SCP7FfK/35zAlDP7c/fLa5h856us2KAOFCLHgx3uJrCZ/b27P3uIdf/o7n8IrbJjVF5e7hUVFVGXkVTuzrifzGH0gG78+vNnRF1Oh/eXFZv418eXsKO+kT7d8vj9NWM4oaRL1GWJpDwzW+Tu5e3bj3QGNdvMXjSz9r3vAG5JTmlyKO/W7mbDjr1MGKL7T3Fw4Um9ePamczhjYBE1Oxv49F2vMnd1+vYeFYnakQLqbeBhYL6ZXd5unQaEC9krwZef7j/FR48uuTwy9Wyeuekc+hR24su/W8h9r6k7ukgYjhRQ7u6/BS4E/tXMfmdmnfetC7c0eaVyMwO6d2ZAj85H3liOGzOjrLgLj39jHJ8Y1pNbZy3je08upbG5NerSRNLKUXWScPdVwFhgE/CmmZ0ValVCS6uz4L2tjB/S/u+WJS665GYx7Ytn8PXzB/PI6+/z+Xvms2nn3qjLEkkbRwqo/Zfx3L3Z3W8m8UezjwBDwyyso1u5cRe79jZzpqZ2j7WMDONfJw7nF1NGsbR6J5/85Vxeq9TfTIkkw5EC6kftG9z9JeAM4MdhFCQJC9duBVBApYjJo/ox64bxFHbK5gv3LuDOOatpbdVVcJGP47AB5e5PHaJ9m7v/JJSKBEgEVJ/CPEqLOkVdihylob26MuuGCVx2al/+67lV/NN9C9lc96EhJEXkKMV6LL6Oyt1ZuHYrZ5Z11+y5KSY/N4tfTBnFf0w+mdcqtzDx53N58Z2aqMsSSUkKqBhav3UPm3Y2cKYmJ0xJZsYXx5Yx65vj6ZGfwzW/X8gPnl6qeaZEjpECKoZeD+4/jdH9p5Q2vHcBT98wnmsnDOL+eeu47FevsLR6R9RliaQMBVQMLXxvK4WdshnaU8PopLq87Ey+f9kIHrh2DLv2NjH5rle5/Zl3dDYlchQUUDGUuP9UREaG7j+li3OGlvDcTedx+eml/Oald7n0F3N5/b2tUZclEmsKqJjZUtfAms27KdflvbRT2Dmb2y8/lQevPYum1lY+d/c8/v2pJeza2xR1aSKxpICKmberEvcoRvfvFm0hEpoJQ4t59qZzuXbCIB5e8D4X/vdfeerNao3nJ9KOAipmFq/fToah+Z/SXOecLL5/2Qie/MZ4+hTmcdOji/nc3fNY9oE6UYjso4CKmbeqtnNir67k52ZFXYocB6f178aT3xjP7Z89hXdrd/OpX73C959ayvb6xqhLE4mcAipG3J231m/ntNJuUZcix1FGhnHFmQN48f+czxfPHshDC9Zx7k9f5O6/vqveftKhKaBiZP3WPWyrb+I03X/qkAo7Z/OjySOZfWNiUsT//PM7XPBfL/GHivW0aFw/6YBCDSgzm2hmK82s0sxuPsj6XDN7NFi/wMzKgvYewUy+dWZ2Z5g1xsniqu0AnNZf9586suG9C/jdNWN45KtnU9I1l3+Z+TaX/OJlnlu2UR0ppEMJLaDMLBO4C7gEGAFcaWYj2m12LbDN3YcAdwC3B+17ge8D/xxWfXH01vrt5GVncGKvrlGXIjEwdnAPnrp+PHdddTqNza1MfWARn/zlKzyzdKNGSpcOIcwzqDFApbuvcfdGYAYwud02k4H7gtczgQvNzNx9t7u/QiKoOoy31m9nZN9CsjN15VUSzIxPntqHF75zHv/1j6dR39jM1x5cxKW/nMvsJRsUVJLWwvwm7Aesb7NcFbQddBt3bwZ2AEc9hayZTTWzCjOrqK2t/ZjlRquppZWlH+zQ/Sc5qKzMDC4/o5QXvnMed1xxGo0trXzjoTe46I6/8sjr76szhaSllP6p7u7T3L3c3ctLSkqiLudjWblxF3ubWhVQclhZmRn8w+hSnv/2efzyytF0ysnklieWMP4nc/j5C6vYovmnJI2E+cc21UD/NsulQdvBtqkysyygENgSYk2xte8PNE/RH+jKUcjMMCad1pdPndqHeWu2cM/c9/j5C6v5zUvv8tkzSvnS2IEM710QdZkiH0uYAbUQGGpmg0gE0RTgqnbbzAKuBuYBlwNzvIN2U1r+wU7yczIZ2L1z1KVICjEzxg0uZtzgYiprdnHvK+8xc1EVDy94n/KBRXz+7AFcMrIPedmZUZcqcswszDwws0uBnwOZwHR3/7GZ3QZUuPssM8sDHgBGA1uBKe6+Jth3LVAA5ADbgYvdffmhPqu8vNwrKipCO5aw/eP/vIY7zPz6uKhLkRS3bXcjMxdV8dCCdazdUk9R52z+sbw/U87szwklmsJF4sfMFrl7+Yfa0+WEJZUDqrXVOfVHz/GZ0/tx2+SRUZcjaaK11Zm3ZgsPzl/Hc8s30dLqnD6gG585vZTLTu1Dt845UZcoAhw6oDTgWwys31ZPXUMzI/ronoEkT0aGMX5IMeOHFFOzcy9PvlnN429U8e9PLeW2Py7nwpN68pnTSznvxBJyslK6v5SkKQVUDCz7YCcAI/oqoCQcPQvyuO68wUw99wSWfbCTJ96oZtZb1fx56Ua6dc7m4hG9uPSUPowbXKywkthQQMXA8g92kplhGkFCQmdmjOxXyMh+hdxy6XDmrq7lj29t4M9LNvJYRRUFeVn83YheXDqyDxOGFqtzhURKARUDyzfsZEhJF30ZyHGVnZnBBcN7ccHwXjQ0t/Bq5Wb+9PZGnl++kSfeqCY/J5NzhpZwwfCenD+8hJ5d86IuWToYBVQMLP9gJ2MHH/UAGiJJl5uVuT+sGptP4bV3N/Pssk28+E4NzyzbCCT+Ru8Tw3tywfCenNqvkIwMi7hqSXcKqIhtqWtg48696iAhsZGTlcH5w3py/rCeuDsrNuzixZU1zHmnhjvnrOaXf1lNUedsxg7uwdjBxYwf3INBxfmYKbAkuRRQEVuxYRegDhIST2bGiL4FjOhbwPWfGMK23Y38dVUtc1dv5rV3NzN7SeLsqndBHuMG92DckGLGlHWnf/dOCiz52BRQEVu+ITHEkc6gJBUU5efw6dH9+PTofrg767bU8+q7m3nt3S28tKqWJ95MjGZW0jWXMwYUUV5WxOkDixjZt1C9A+WYKaAi9s7GXfQqyKUoX380KanFzCgrzqesOJ/PnzWQ1lZnVc0uKtZuY9G6xGPf/aucrAxOKy1k9IAiRvYr5JR+hQzs3ln3seSwFFARW72pTt3LJS1kZBjDexcwvHcBXzh7IAA1O/fyxvvbqFi7jYp12/j9a2tpbG4FoGtuFif3K+CUfoWcUtpNoSUfooCKUGurU1lTx5VjBkRdikgoehbkMXFkHyaO7AMk5j1btWkXS6p2sKR6B0urd3Dfa+tobHkPgM45mQzt1ZVhvbowrHcBw3p1ZVjvrhR3ydE9rQ5IARWhqm172NPUwom9NICndAzZmRmc3LeQk/sWMiVoa2xOhNbS6h28s3EXqzbt4i8ranisomr/ft3zcxjWqytDe3XhhOJ8BpUknvt260SmzrjSlgIqQqs2JXrwDdUlPunAcrIy9o9u0dbmugZWbdy1P7Te2biLJ96opq6h+W/7ZmYwsEdnBhXnM6gknxOK8xnYI5/Sok70LsgjK1MdM1KZAipCK4OA0hmUyIcVd8mleEgu44YU729zdzbXNfLe5t28t7mONZt3817tbt7bvJuXVtbS2NK6f9usDKNPtzxKu3WmtKgTpUX7njtR2r0zvbrmKsBiTgEVodWbdtG3MI+uedlRlyKSEsyMkq65lHTNZcyg7gesa2l1qrft4f2t9VRtq6dq2x6qttWzftseXl5dy6adDQdsn2GJEOxdmEevgjx6FeTSuyCPngV59C5ItPUuyKOgU5buf0VEARWhVZvqdHlPJEkyM4wBPTozoMfBZ6VuaG7hg+1794fXB9v3sGnnXjbubGD91noWrt3K9vqmD+2Xl51Bj/xcirvk0KNLLt3zc+jRJYfi/Davg/bu+TkaUzOJFFARaWl13q2tY/wQjcEncjzkZmUm7lUV5x9ym71NLdTsbGDTrr1s3LGXTTsTjy11jWze3UjNrr2s2LCTLXWNB1xObKtLbhZF+dkUdko8unXKoaBTm+XOf3u9/9E5m665OlNrTwEVkfe31tPQ3KozKJEYycvOPOxZ2D7uTl1DM1vqGtmyuyF4bmRLXQOb6xrZsaeJ7fWJ5407drJjTxM79jTR1HLoGcwzDPJzsuiSl0V+bhZd2jwSy5n713XNPXCb/NwsOuVk0ik7k7zsTDrnJJ5TvYdjqAFlZhOBXwCZwD3u/pN263OB+4EzgC3AFe6+Nlh3C3At0AJ8y92fDbPW423V/g4SCiiRVGNmdM3LpmteNmWHOSNry93Z09TC9vqm/YG1vb6JnXv+tlzX0ExdQzO7g+e6hmZqdu2lbm/Q3thCS+uhQ669nKwMOmUngqtTEFqdsjMOCLO26/KyMsjJyiA3KzN4/tvy314fuE3XvCx6FYQzFUtoAWVmmcBdwEVAFbDQzGa5+/I2m10LbHP3IWY2BbgduMLMRgBTgJOBvsALZnaiu7eEVe/xtnpfF/Oe6sEn0hGYGZ1zsuick0Xfbp0+0nu4O3ubWg8Isl17E897m1vY09jC3qYW9jS1sKexlT1NwXJj0NZmeXt9U2I5WFff2EJjSyt+9PkHwKj+3Xjq+vEf6XiOJMwzqDFApbuvATCzGcBkoG1ATQZ+GLyeCdxpiYuwk4EZ7t4AvGdmlcH7zQur2B/9cRnLg6nXj4faXQ2UFnUiP1dXWUXk6JhZ4uwnJ5OSrrlJf393p6nFaWxppaGphYbmVhqbW9s8t+xfbmhOrC8IsRdymN+O/YD1bZargLMOtY27N5vZDqBH0D6/3b792n+AmU0FpgIMGJBawwVlZhhn9u9+5A1FRI4TMyMny8jJyqBLDH48R1/Bx+Du04BpAOXl5cd4YnqgWz91clJqEhGR5Ajzz6irgf5tlkuDtoNuY2ZZQCGJzhJHs6+IiKSxMANqITDUzAaZWQ6JTg+z2m0zC7g6eH05MMfdPWifYma5ZjYIGAq8HmKtIiISM6Fd4gvuKd0APEuim/l0d19mZrcBFe4+C7gXeCDoBLGVRIgRbPcYiQ4VzcD16dSDT0REjsz8WPsUxpSZ1QLrPubbFAObk1BOXOn4Ulc6Hxvo+FLdxz2+ge5e0r4xbQIqGcyswt3Lo64jLDq+1JXOxwY6vlQX1vFprHkREYklBZSIiMSSAupA06IuIGQ6vtSVzscGOr5UF8rx6R6UiIjEks6gREQklhRQIiISSwqogJlNNLOVZlZpZjdHXU8ymdl0M6sxs6VR15JsZtbfzF40s+VmtszMboy6pmQyszwze93M3gqO70dR15RsZpZpZm+a2f9GXUuymdlaM1tiZovNrCLqepLNzLqZ2Uwze8fMVpjZ2KS+v+5B7Z+7ahVt5q4Crmw3d1XKMrNzgTrgfncfGXU9yWRmfYA+7v6GmXUFFgGfTqN/OwPy3b3OzLKBV4Ab3X3+EXZNGWb2HaAcKHD3y6KuJ5nMbC1Q7u5p+Ue6ZnYfMNfd7wmGtOvs7tuT9f46g0rYP3eVuzcC++auSgvu/jKJoaTSjrtvcPc3gte7gBUcZGqWVOUJdcFidvBIm1+VZlYKfBK4J+pa5NiYWSFwLokh63D3xmSGEyig9jnY3FVp8yXXUZhZGTAaWBBxKUkVXAJbDNQAz7t7Oh3fz4HvAq0R1xEWB54zs0XB/HXpZBBQC/wuuER7j5nlJ/MDFFCSFsysC/A4cJO7H7+pkY8Dd29x91Ekpp0ZY2ZpcZnWzC4Datx9UdS1hGiCu58OXAJcH1xuTxdZwOnAb9x9NLAbSOr9ewVUguafSmHBvZnHgYfc/Ymo6wlLcPnkRWBixKUky3hgUnCfZgZwgZk9GG1JyeXu1cFzDfAkidsJ6aIKqGpzRj+TRGAljQIq4WjmrpIYCjoR3AuscPefRV1PsplZiZl1C153ItGR551Ii0oSd7/F3UvdvYzE/3Nz3P0LEZeVNGaWH3TcIbj0dTGQNj1p3X0jsN7MhgVNF5KYIilpUnrK92Q51NxVEZeVNGb2CHA+UGxmVcCt7n5vtFUlzXjgi8CS4D4NwPfcfXZ0JSVVH+C+oKdpBvCYu6ddd+w01Qt4MvEbiizgYXd/JtqSku6bwEPBD/s1wDXJfHN1MxcRkVjSJT4REYklBZSIiMSSAkpERGJJASUiIrGkgBIRkVhSQIlEwMx6BCNcLzazjWZWHbyuM7NfR12fSByom7lIxMzsh0Cdu/9X1LWIxInOoERixMzO3zcvkpn90MzuM7O5ZrbOzD5jZj8N5hd6JhjiCTM7w8z+GgxI+mwwBYlIylNAicTbYOACYBLwIPCiu58C7AE+GYTUr4DL3f0MYDrw46iKFUkmDXUkEm9/dvcmM1tCYhiufUPlLAHKgGHASOD5YEidTGBDBHWKJJ0CSiTeGgDcvdXMmvxvN41bSfz/a8Ayd0/qVNsicaBLfCKpbSVQYmZjITH1iJmdHHFNIkmhgBJJYe7eCFwO3G5mbwGLgXGRFiWSJOpmLiIisaQzKBERiSUFlIiIxJICSkREYkkBJSIisaSAEhGRWFJAiYhILCmgREQklv4/V2ufcbMkEaQAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "plt.subplot( 2, 1, 1)\n", "plt.plot( [sol1.t[0],time_x_on,time_x_on,time_x_off,time_x_off,sol3.t[-1]], [0,0,1,1,0,0] )\n", "plt.ylabel( \"X\")\n", "plt.subplot( 2, 1, 2 )\n", "plt.plot(sol1.t, sol1.y[1,:].T); #Simulation before signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol2.t, sol2.y[1,:].T); #Simulation after signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol3.t, sol3.y[1,:].T); #Simulation after signal on\n", "plt.xlabel( 'Time')\n", "plt.ylabel( \"Z\")\n", "plt.tight_layout()" ] }, { "cell_type": "code", "execution_count": 40, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa8AAAEYCAYAAADrpHnMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAAAs70lEQVR4nO3deZgV5Zn///fdO91A09DN2kAD4gKiqC2IoDHGBY2CoyZRv2ZIxkRnRqIZM5NoZiYac81vssxk1TgyqGPcUFEjRlyIYlQ2aRBBUKBBkEaWZqeB3u/fH6cwx7ZZlFNU1+nP67rOdU49VXXOXQJ+Tj31nHrM3REREYmTjKgLEBER+awUXiIiEjsKLxERiR2Fl4iIxI7CS0REYicr6gJSpbi42MvKyqIuQ0REUmTBggVb3L2ktXVpE15lZWVUVFREXYaIiKSIma090LrQug3N7H4z22xm7x5gvZnZb82s0swWm9mpSesmmNnK4DEhrBpFRCSewrzm9X/A2IOsvwgYHDyuB+4BMLOuwO3ASGAEcLuZFYVYp4iIxExo3Ybu/rqZlR1kk/HAHzxxi4+5ZtbFzHoB5wAz3H0bgJnNIBGCj4VVa1Qenfchzy5aH3UZcoTGD+/DNSP7RV2GSLsS5WjDPsC6pOWqoO1A7Z9iZtebWYWZVVRXV4dWaFieXbSeZRt2RV2GHIFlG3bpC4hIBGI9YMPdJwGTAMrLy2N5k8YhvTrz+A2joi5DPqev3Tsn6hJE2qUoz7zWA32TlkuDtgO1i4iIANGG1zTgb4NRh2cAO919A/AScIGZFQUDNS4I2kRERIAQuw3N7DESgy+KzayKxAjCbAB3/x9gOnAxUAnsBb4ZrNtmZj8B5gdvdef+wRsiIiIQ7mjDqw+x3oEbD7DufuD+MOoSEZH4070NRUQkdhReIiISOwovERGJHYWXiIjEjsJLRERiR+ElIiKxo/ASEZHYUXiJiEjsKLxERCR2FF4iIhI7Ci8REYkdhZeIiMSOwktERGJH4SUiIrGj8BIRkdgJNbzMbKyZLTezSjO7tZX1vzKzRcFjhZntSFrXlLRuWph1iohIvIQ5k3ImcDdwPlAFzDezae6+bP827v5PSdt/Bzgl6S32ufvwsOoTEZH4CvPMawRQ6e6r3b0emAKMP8j2VwOPhViPiIikiTDDqw+wLmm5Kmj7FDPrDwwAXk1qzjOzCjOba2aXHWC/64NtKqqrq1NUtoiItHVtZcDGVcBUd29Kauvv7uXANcCvzWxQy53cfZK7l7t7eUlJydGqVUREIhZmeK0H+iYtlwZtrbmKFl2G7r4+eF4NvMYnr4eJiEg7FmZ4zQcGm9kAM8shEVCfGjVoZscDRcCcpLYiM8sNXhcDo4FlLfcVEZH2KbTRhu7eaGYTgZeATOB+d19qZncCFe6+P8iuAqa4uyftfgJwr5k1kwjYnyaPUhQRkfYttPACcPfpwPQWbT9qsXxHK/vNBoaFWZuIiMRXWxmwISIictgUXiIiEjsKLxERiR2Fl4iIxI7CS0REYkfhJSIisaPwEhGR2FF4iYhI7Ci8REQkdhReIiISOwovERGJHYWXiIjEjsJLRERiR+ElIiKxE2p4mdlYM1tuZpVmdmsr679hZtVmtih4fCtp3QQzWxk8JoRZp4iIxEto83mZWSZwN3A+UAXMN7NprUwq+bi7T2yxb1fgdqAccGBBsO/2sOoVEZH4OOCZl5n1O8i6sw7jvUcAle6+2t3rgSnA+MOs60JghrtvCwJrBjD2MPcVEZE0d7Buw9fM7PvBGRQAZtbDzB4GfnUY790HWJe0XBW0tXSFmS02s6lm1vez7Gtm15tZhZlVVFdXH0ZJIiKSDg4WXqcBg4BFZnaumd0MvAXMIXFWlQrPAWXufhKJs6sHP8vO7j7J3cvdvbykpCRFJYmISFt3wGteQXfdDUFo/Rn4CDjD3asO873XA32TlkuDtuTP2Jq0OBn4edK+57TY97XD/FwREUlzB7vm1cXM7gW+SeJ601TgBTM79zDfez4w2MwGmFkOcBUwrcVn9EpaHAe8F7x+CbjAzIrMrAi4IGgTERE56GjDhcDvgRvdvRF42cyGA783s7XufvXB3tjdG81sIonQyQTud/elZnYnUOHu04CbzGwc0AhsA74R7LvNzH5CIgAB7nT3bZ/7KEVEJK0cLLzObtlF6O6LgDPN7NuH8+buPh2Y3qLtR0mvbwNuO8C+9wP3H87niIhI+3LAbsODXdty9/8NpxwREZFD0+2hREQkdhReIiISOwovERGJHYWXiIjEjsJLRERiR+ElIiKxo/ASEZHYUXiJiEjsKLxERCR2FF4iIhI7Ci8REYkdhZeIiMSOwktERGJH4SUiIrETaniZ2VgzW25mlWZ2ayvrbzGzZWa22MxeMbP+SeuazGxR8JjWcl8REWm/DjYZ5RExs0zgbuB8oAqYb2bT3H1Z0mZvA+XuvtfM/gH4OfC1YN0+dx8eVn0iIhJfYZ55jQAq3X21u9cDU4DxyRu4+0x33xsszgVKQ6xHRETSRJjh1QdYl7RcFbQdyHXAC0nLeWZWYWZzzeyy1nYws+uDbSqqq6uPuGAREYmH0LoNPwszuxYoB76Q1Nzf3deb2UDgVTNb4u6rkvdz90nAJIDy8nI/agWLiEikwjzzWg/0TVouDdo+wczOA/4VGOfudfvb3X198LwaeA04JcRaRUQkRsIMr/nAYDMbYGY5wFXAJ0YNmtkpwL0kgmtzUnuRmeUGr4uB0UDyQA8REWnHQus2dPdGM5sIvARkAve7+1IzuxOocPdpwC+AjsCTZgbwobuPA04A7jWzZhIB+9MWoxRFRKQdC/Wal7tPB6a3aPtR0uvzDrDfbGBYmLWJiEh86Q4bIiISOwovERGJHYWXiIjEjsJLRERiR+ElIiKxo/ASEZHYUXiJiEjsKLxERCR2FF4iIhI7Ci8REYkdhZeIiMSOwktERGJH4SUiIrGj8BIRkdgJNbzMbKyZLTezSjO7tZX1uWb2eLB+npmVJa27LWhfbmYXhlmniIjES2jhZWaZwN3ARcAQ4GozG9Jis+uA7e5+DPAr4GfBvkNIzLw8FBgL/D54PxERkVAnoxwBVLr7agAzmwKMB5JnRB4P3BG8ngrcZYkplccDU9y9DvjAzCqD95sTVrE/fm4pyz7aFdbbt2rZhl0M6dX5qH6mpN6yDbv42r2h/dUUiaUhvTtz+6VDQ3v/MMOrD7AuabkKGHmgbdy90cx2At2C9rkt9u3T8gPM7HrgeoB+/fqlrPCjZUivzowf/qnDkhjRn59INMIMr9C5+yRgEkB5ebkfyXuF+Q1B0tc1I/txzcj4fXESibswB2ysB/omLZcGba1uY2ZZQCGw9TD3FRGRdirM8JoPDDazAWaWQ2IAxrQW20wDJgSvrwRedXcP2q8KRiMOAAYDb4VYq4iIxEho3YbBNayJwEtAJnC/uy81szuBCnefBtwHPBQMyNhGIuAItnuCxOCORuBGd28Kq1YREYkXS5zoxJ+ZVQNrj/BtioEtKSinLdMxxl+6Hx/oGNPFkR5jf3cvaW1F2oRXKphZhbuXR11HmHSM8Zfuxwc6xnQR5jHq9lAiIhI7Ci8REYkdhdcnTYq6gKNAxxh/6X58oGNMF6Edo655iYhI7OjMS0REYkfhJSIisaPwChxq7rG4M7P7zWyzmb0bdS1hMLO+ZjbTzJaZ2VIzuznqmlLNzPLM7C0zeyc4xh9HXVMYzCzTzN42sz9FXUsYzGyNmS0xs0VmVhF1PWEwsy5mNtXM3jez98xsVMo/Q9e8Pp57bAVwPok72M8Hrnb3ZQfdMUbM7GygBviDu58YdT2pZma9gF7uvtDMOgELgMvS7M/QgAJ3rzGzbOBN4GZ3n3uIXWPFzG4ByoHO7n5J1PWkmpmtAcrdPW1/oGxmDwJvuPvk4PaA+e6+I5WfoTOvhI/nHnP3emD/3GNpw91fJ3ELrrTk7hvcfWHwejfwHq1MoxNnnlATLGYHj7T69mlmpcCXgclR1yKfj5kVAmeTuP0f7l6f6uAChdd+rc09llb/42tPzKwMOAWYF3EpKRd0qS0CNgMz3D3djvHXwPeB5ojrCJMDL5vZgmBOwnQzAKgGHgi6fyebWUGqP0ThJWnFzDoCTwHfdfejOzX2UeDuTe4+nMQ0QSPMLG26gM3sEmCzuy+IupaQjXH3U4GLgBuDLv10kgWcCtzj7qcAe4CUjyNQeCVo/rA0EFwHegp4xN2fjrqeMAXdMDOBsRGXkkqjgXHBNaEpwLlm9nC0JaWeu68PnjcDz5C4bJFOqoCqpF6BqSTCLKUUXgmHM/eYtGHBYIb7gPfc/ZdR1xMGMysxsy7B6w4kBhi9H2lRKeTut7l7qbuXkfg3+Kq7XxtxWSllZgXBgCKCrrQLgLQaAezuG4F1ZnZc0PQlEtNbpVRo83nFyYHmHou4rJQys8eAc4BiM6sCbnf3+6KtKqVGA18HlgTXhAB+6O7Toysp5XoBDwajYzOAJ9w9LYeTp7EewDOJ71pkAY+6+4vRlhSK7wCPBCcDq4FvpvoDNFReRERiR92GIiISOwovERGJHYWXiIjEjsJLRERiR+ElIiKxo/ASEZHYUXiJiEjsKLxERCR2FF4iIhI7Ci8REYkdhZeIiMRO2tyYt7i42MvKyqIuQ0REUmTBggVb3L2ktXVpE15lZWVUVFREXYaIiKSIma090Dp1G4qISOykzZmXRMfdcQczCOYpEhEJlcJLDpu7U7m5hjmrt/LOup0s37SLTbvq2FpTR3MwLVxBTiYlnXLpU9SB43p0ZmjvzpwxqBt9unSItngRSSsKLzmk6t11PDJvLdPe+YjV1XsAKO6Yw5DehZzYu5DijrlkZRrusLu2kc27a1m3bS+PvrWW2oZmAAaWFHD24BLGDe/NKX276AxNRI6IwksOaMPOffx6xkqeeXs9Dc3NjBrYjW+eWcYXj+9Ony4dDhlATc3Oys27eXPlFt6s3MJjb33I/81eQ/9u+fzNKX24ZmQ/unfKO0pHIyLpxNw96hpSory83DXaMDVqG5r439dX8/vXVtHkzlfLS/m70QMYWNLxiN53d20DL767kT8uWs+syq1kZxqXnNSb68YM4MQ+hSmqXkTShZktcPfyVtcpvCRZ5ebdTHz0bd7fuJuLh/XktotOoG/X/JR/zgdb9vDg7DVMXVBFTV0j5w/pwXfPG8zQ3goxEUlQeMlhebJiHf/+7LsU5GTxi6+cxLnH9wj9M3fVNvDAm2uY/OZqdtc2csGQHvzgouMZdIRneSISfwovOSh351czVvDbVys5c1A3fv214XTvfHSvRe3c18ADsz5g8hsfUNvQxIQzy7jp3MEU5mcf1TpEpO1QeMkBNTU7P3hqMVMXVPHV8lL+42+GkZ0Z3W/Xq3fX8csZy5kyfx1dOmRzy/nHcs3I/mRmaHSiSHtzsPDSHTbaseak4PrueYP52RUnRRpcACWdcvnPy0/iT98Zw3E9O/Hvzy7l8ntm896GXZHWJSJtS6j/pzKzsWa23MwqzezWVtbfYmbLzGyxmb1iZv2T1k0ws5XBY0KYdbZH7s6Ppr37cXB997xj29Rvr4b2LuSxb5/Bb64aTtW2vVz6uzf52YvvU9vQFHVpItIGhBZeZpYJ3A1cBAwBrjazIS02exsod/eTgKnAz4N9uwK3AyOBEcDtZlYUVq3t0d0zK3l47ofc8IWB3PylwVGX0yozY/zwPrzyvS9w+al9uOe1VVzwq9eZXbkl6tJEJGJhnnmNACrdfbW71wNTgPHJG7j7THffGyzOBUqD1xcCM9x9m7tvB2YAY0OstV158d2N/NfLK7hseG9uHXt8mzrjak2X/Bx+fuXJPPrtkWRmGNdMnsePn1uqszCRdizM8OoDrEtargraDuQ64IXPsq+ZXW9mFWZWUV1dfYTltg/vb9zFLU8s4uTSQn56xUltPriSnTmomOk3ncU3zizjgVlr+PJv3+CddTuiLktEItAmBmyY2bVAOfCLz7Kfu09y93J3Ly8paXW+Mkmyr76JGx9ZSEFuFpP+tpy87MyoS/rMOuRkcse4oTzyrZHsrW/i8ntm88sZK2hoao66NBE5isIMr/VA36Tl0qDtE8zsPOBfgXHuXvdZ9pXP5ifPL2NV9R5+9dXh9DjKv+NKtdHHFPPid89m/PDe/PaVlVx5z2zWbt0TdVkicpSEGV7zgcFmNsDMcoCrgGnJG5jZKcC9JIJrc9Kql4ALzKwoGKhxQdAmn9PLSzfy6LwPueHsgYwZXBx1OSlR2CGbX351OPf8v1P5YMsevvzbN3l2kb7jiLQHoYWXuzcCE0mEznvAE+6+1MzuNLNxwWa/ADoCT5rZIjObFuy7DfgJiQCcD9wZtMnnsGNvPbc9vYShvTvzvQuOi7qclLtoWC+m33wWx/XsxM1TFvH9qe+wt74x6rJEJES6w0Y78IOpi5m6sIrnJo5hSO/OUZcTmsamZn7955Xc/VolA4sLuOuaUzmhV/oer0i60x022rF5q7fyeMU6vjVmQFoHF0BWZgb/fOFxPHLdSHbVNjL+7lk8NGcN6fIFTUT+SuGVxuobm/nhM0soLerAzee1zR8ih+HMY4p54eazOHNQN/792aX8w8ML2bmvIeqyRCSFFF5p7JF5a1lVvYcfjxtKfk77mjS7uGMu9084nR9efDx/fm8Tl/7uTZZU7Yy6LBFJEYVXmtq5t4HfvLKS0cd049zju0ddTiQyMozrzx7E4zeMorGpmSvumc0f1I0okhYUXmnqrpkr2bmvgX+9eEis7qIRhtP6F/H8TWcx+phu/OjZpUx89G121aobUSTOFF5paN22vTw4ey1Xnlqa9oM0DldRQQ73TTidWy86nheXbuTS373Ju+vVjSgSVwqvNPS7V1eCkZa/6ToSGRnG339hEFOuP4O6hmYuv2c2D89dq25EkRhSeKWZD7fu5emF67lmRD96Fsb7FlBhOb2sK9NvPotRA7vxb398l+889ja71Y0oEisKrzRz98xKMjKMfzhnUNSltGldC3J44Bun8/2xx/HCuxsZd9csln2k2ZpF4kLhlUbWbdvLUwuruGZEv9jfePdoyMgw/vGcY3js22ewt76Ry34/i0fnfahuRJEYUHilkXv+surj6zpy+EYM6MrzN53FyAFd+eEzS/ju44vYU6d7I4q0ZQqvNLG1po6nFlRxxal9dK3rcyjumMuD3xzBP19wLM+98xGX3vUm729UN6JIWxVqeJnZWDNbbmaVZnZrK+vPNrOFZtZoZle2WNcU3Gn+47vNy4E9PPdD6hqbuW7MwKhLia2MDGPiuYN55FtnsLu2kfF3zeLx+epGFGmLQgsvM8sE7gYuAoYAV5vZkBabfQh8A3i0lbfY5+7Dg8e4VtZLoLahiYfmruHc47tzTPeOUZcTe6MGdWP6TWdxellXfvDUEr73xDvqRhRpY8I88xoBVLr7anevB6YA45M3cPc17r4Y0BzuR+CPb69nS0093xozIOpS0kZJp1we/LsR/NN5x/LMovWMu+tNlm/cHXVZIhIIM7z6AOuSlquCtsOVZ2YVZjbXzC5rbQMzuz7YpqK6uvoISo0vd+f+WR8wpFdnRg3qFnU5aSUzw7j5vME8ct1Idu5rZPzdb/JkxbpD7ygioWvLAzb6B5OQXQP82sw+NYTO3Se5e7m7l5eUlBz9CtuA+Wu2s2JTDd84s6zd38MwLGceU8z0m8dwSt8i/mXqYr73hGZqFolamOG1HuibtFwatB0Wd18fPK8GXgNOSWVx6eKReWvplJfFJSf3irqUtNa9Ux4Pf2skN31pME+/XcX4u2axcpO6EUWiEmZ4zQcGm9kAM8sBrgIOa9SgmRWZWW7wuhgYDSwLrdKY2lpTxwtLNnLFqaXtbr6uKGRmGLecfywP/d1Itu+tZ9xds3jsrbV8sGVP1KWJtDuhhZe7NwITgZeA94An3H2pmd1pZuMAzOx0M6sCvgLca2ZLg91PACrM7B1gJvBTd1d4tTB1QRX1Tc1cM7Jf1KW0K2MGF/P8TWdxUmkhtz39Lpf87g39JkzkKLN0+Q1LeXm5V1RURF3GUdPc7Hzxv1+jR6c8nvj7UVGX0y41NjXzyxkreGphFfWNzTx03UhO7FMYdVkiacPMFgRjHz6lLQ/YkIOY98E21m7dy9Uj+x56YwlFVmYG3x97PE/cMIr8nCyu+d+5LFq3I+qyRNoFhVdMPb2wio65WYwdqoEaUevfrYDHbziDLvk5XDt5HgvWbou6JJG0p/CKob31jUxfsoGLh/WkQ05m1OUIUFqUz+M3nEFJp1y+ft9bzK7cEnVJImlN4RVDLy/dxJ76Ji4/tTTqUiRJr8IOPH79GZQWdeAbD8znT4s/irokkbSl8IqhpxZWUVrUgRFlXaMuRVro3jmPJ284k5P7FvKdx97m/2Z9EHVJImlJ4RUzG3fWMqtyC5ef0oeMDN1Roy0qzM/moetGcv4JPbjjuWX87MX3dWd6kRRTeMXMc+98RLPD36jLsE3Ly87knmtP45qR/bjntVV878l3qG/U/adFUuWg4WVmFx5k3VdSX44cyp8Wf8SwPoUMKC6IuhQ5hMwM4z8uO5Fbzj+Wpxeu59rJ89haUxd1WSJp4VBnXtPNbKaZtXY3+NvCKEgObN22vbxTtZNLTtLw+LgwM2760mB+e/UpvFO1g3F3zdLdOERS4FDhtZjERJFzW850DOiCy1H2p8UbALh4mMIrbsad3JsnbhhFQ1MzV/x+Nn9etinqkkRi7VDh5e7+v8CXgB+Y2QNmlr9/XbilSUvPL/mI4X270Ldr/qE3ljbn5L5dmDZxDIO6d+TbD1Xwu1dW0tysf0Yin8dhDdhw9xXAKGAT8LaZjTyc/cxsrJktN7NKM7u1lfVnm9lCM2tseWZnZhPMbGXwmHA4n5fO1mzZw7vrd6nLMOZ6FubxxA2jGH9yb/57xgr+7sH5bN9TH3VZIrFzqPD6uGvQ3Rvd/VbgBuAxYPBBdzTLBO4GLgKGAFeb2ZAWm30IfINE12Tyvl2B24GRwAjgdjMrOtTBpLPnl6jLMF3kZWfyq68N5yeXncjsyq1c8rs3dU9Ekc/oUOH145YN7v4acBrwH4fYdwRQ6e6r3b0emAKMb/Fea9x9MdByDPGFwAx33+bu24EZwNhDfF5ae/HdjZzSrwu9u3SIuhRJATPj62f058lgRoCv/M9sHpy9Rr8HEzlMBw0vd//jAdq3u/tPD/HefYB1SctVQdvhOJJ9086GnftYsn4nFwzpGXUpkmIn9+3C8zeNYcwxxdw+bSnf/kOFhtOLHIZY/0jZzK43swozq6iuro66nND8+b3NAJw/pHvElUgYuuTncN+E0/n3S4bw+ootXPjrN5i5fHPUZYm0aWGG13ogebKp0qAtZfu6+yR3L3f38pKSks9daFs3Y9kmBhQXMKikY9SlSEgyMozrxgzg2Ymj6VaQwzcfmM/tz75LbUNT1KWJtElhhtd8YLCZDTCzHOAqYNph7vsScIGZFQUDNS4I2tqd3bUNzFm1hfOH9MBMP61Ldyf06syzE0fzzdFlPDhnLRf/5g3mr9H8YCIthRZe7t4ITCQROu8BT7j7UjO708zGAZjZ6WZWBXwFuNfMlgb7bgN+QiIA5wN3Bm3tzusrttDQ5Jx3Qo+oS5GjJC87k9svHcoj3xpJfVMzX713Drc/+y41dY1RlybSZli6jG4qLy/3ioqKqMtIue9OeZu/rKim4t/OJ1N3kW939tQ18ouXlvPgnDX0LuzAf14+jLOPTd8ucpFkZrbA3ctbWxfrARvprqGpmVff38y5x/dQcLVTBblZ3DFuKE/eMIrc7Az+9v63mPjoQjburI26NJFIKbzasPlrtrGrtpHzh6jLsL0rL+vK9JvO4p/OO5YZyzZx7n+/xr1/WaVpVqTdUni1Ya++t5mczAzOPrY46lKkDcjLzuTm8wbz51u+wJmDivnPF97not+8zqzKLVGXJnLUKbzasNdXVjNiQFfyc7KiLkXakL5d85k8oZz7JpTT0OT8v8nz+NaDFVRu3h11aSJHjcKrjdqwcx8rNtXorEsO6Esn9ODlfzqbf7nwOOat3soFv3qd255ewuZduh4m6U/h1Ua9sSLRFaSRZXIwedmZ3PjFY3jtX87hb0eV8WTFOr7wi9f45YwV7K5tiLo8kdAovNqov6yspkfnXI7r0SnqUiQGunXM5Y5xQ/nzLV/g3BO689tXVjLmZzO569WVCjFJSwqvNqip2Xlz5RbOGlyiu2rIZ1JWXMDd15zKcxPHcHpZEf/18grO+vlM7p5ZqRCTtKLwaoMWV+1g574GdRnK5zastJDJE07nuYljKO9fxC9eWs5ZP5/JL2es0F3rJS0ovNqg11dswQzOOkaDNeTI7A+xaRNHU96/K799ZSVn/vRV/vWZJXywZU/U5Yl8bhqD3Qa9vrKak/oUUlSQE3UpkiZOKu3C5AnlVG6uYfIbq3myoopH3/qQC4f05NtnD+TUfl3URS2xovBqY3bua2DRuh384zmDoi5F0tAx3Tvy0ytO4pYLjuUPs9fy0Ny1vLh0Iyf26cy1I/szbnhv/a5QYiHUbkMzG2tmy82s0sxubWV9rpk9HqyfZ2ZlQXuZme0zs0XB43/CrLMtmbNqC03NzlmDdb1LwtO9Ux7/fOFxzL71XH5y2Yk0Njm3Pr2Ekf/fK9wxbSmVm2uiLlHkoEL7imVmmcDdwPlAFTDfzKa5+7Kkza4Dtrv7MWZ2FfAz4GvBulXuPjys+tqqOau20iE7k+F9u0RdirQDBblZfP2M/lw7sh8Va7fz0Jy1PDJvLf83ew0jB3TlitNKuXhYLzrm6mxM2pYw/0aOACrdfTWAmU0BxgPJ4TUeuCN4PRW4y9p5x/vc1dsoLysiJ0tjaeToMTNOL+vK6WVdqd49hCcq1vFkxTq+P3Uxtz+7lLEn9uTK00oZNbAbGZrhQNqAMMOrD7AuabkKGHmgbdy90cx2At2CdQPM7G1gF/Bv7v5Gyw8ws+uB6wH69euX2uojsLWmjuWbdjNueO+oS5F2rKRTLjd+8Rj+8ZxBLPxwO1MXrOdPiz/imbfX07swj8tO6cOXT+rFkF6dNchDItNW+wI2AP3cfauZnQb80cyGuvuu5I3cfRIwCRKTUUZQZ0rNXZ2YLHrUoG6H2FIkfGbGaf27clr/rtx+6RBmLNvEUwuruPf11fz+tVWUdcvn4mG9uHhYL4b2VpDJ0RVmeK0H+iYtlwZtrW1TZWZZQCGw1RPTO9cBuPsCM1sFHAuk31TJSeau3kpBTibD+hRGXYrIJ+RlZ3Lpyb259OTebK2p4+Vlm5i+ZMPHQdY/CLLzTujB8L5dNHmqhC7M8JoPDDazASRC6irgmhbbTAMmAHOAK4FX3d3NrATY5u5NZjYQGAysDrHWNmHO6q2Ul3UlO1PXu6Tt6tYxl6tH9OPqEf3Ytqeel5du5PklG5j0+mrueW0VXQtyOOfYEs49oTtnDS6hsEN21CVLGgotvIJrWBOBl4BM4H53X2pmdwIV7j4NuA94yMwqgW0kAg7gbOBOM2sAmoG/d/dtYdXaFmzeXUvl5hquPK006lJEDlvXghyuGtGPq0b0Y+feBl5fWc2r72/m1eWbefrt9WRlGOVlRXzxuO6MGVzMCT07a8CHpIQleujir7y83Csq4tur+Nw7H/Gdx97m2RtHc7KGyUvMNTU7i9Zt55X3NvPq+5t5f2NiosyuBTmMGtiNM4/pxuhBxfTvlq9rZXJAZrbA3ctbW9dWB2y0O3NWb6VjbhZDe3eOuhSRI5aZ8dfBHt8fezwbd9Yyq3ILs1ZtYXblVp5fsgGAPl06cMbAbpxeVkR5WRGDSjoqzOSwKLzaiLmrtzJiQFeydL1L0lDPwjyuOK2UK04rxd1ZvWUPsyu3MKtyKzOXb+aphVUAdMnPprx/Eaf170p5WRHD+hSSl50ZcfXSFim82oDNu2pZXb2Hq07ve+iNRWLOzBhU0pFBJR35+qiyj8NswZrtzF+zjQVrt/Pn9zYDkJ1pnNCrM8P6FHJSaSHD+nRhcI+OGtQkCq+2YP6a7QCMGKDfd0n7kxxmXw2+wG2pqWPB2u0sXLudJet3Mm3RRzwy70MAcrMyGNK7Myf1KeTEPoWc0Kszx3TvqDO0dkbh1QZUrN1GXnaGrneJBIo75nLh0J5cOLQnAM3Nztpte1lctYMlVTtZvH4nUxdU8eCctQBkWGIW6eN7duK4Hp05rmcnju/ZiX5d8zW6MU0pvNqABWu3c3JpF3WFiBxARoYxoLiAAcUFjB/eB0iMaFyzdQ/LN+7m/Y27Wb5xF0s/2sUL725k/yDqDtmZDOpewIDijgwsLmBgSQEDizsyoKRANxuOOf3pRWxvfSNLP9rF339hYNSliMRKZsZfuxsvHtbr4/a99Y2s2FTD8o27eH/jblZV72HRuu38afFHJP8yqHunXAYUFzCwpCP9u+XTtyif0qIO9O2aT1F+tkY9tnEKr4gtWreDpmanvH/XqEsRSQv5OVkM79vlU9MK1TY08eG2vayurmH1lj2srt7DB1v28OK7G9i+t+ET2xbkZFJalE/frh0oDUKttCifXoV59CzMo7hjrm6BFTGFV8QWBIM1Tu1XFHElIuktLzuTY3t04tgenT61bldtA1Xb9rFu+16qtu9j3bbEc9X2vcxZtZU99U2f2D4zw+jeKZcenfPo2TkRaD0659GzMJeenTvQo3MuxZ1y6ZSbpTO4kCi8IlaxdjvH9uhIYb7u/yYSlc552Qzpnc2QVgZNuTs79jZQtX0fG3buY9OuWjbuqmXjzjo27aqlsrqGWZVb2F3X+Kl9c7IyKC7IoVvHXLp1zKE4eC4JnrsV5H7cVtghWyMmPwOFV4Sam52FH27nkpM0f5dIW2VmFBXkUFSQw7DSA8/4sKeukY27atm0s5ZNu2vZsrueLXvq2FpTz9aaOrbU1LNi42621NRT39Tc6nvkZWfQpUMOXfKzKeyQTZf87L8uJ73u0iGbzh2y6ZSXRcfcLDrmZZGb1b6CT+EVoRWbd7O7tpHy/uoyFIm7gtysjweQHIy7s7uuMSnU6ti6p54dexvYua+BHXsTr3fsa2DNlr3s2LeD7XsbqG9sPfD2y8nMoOP+MAsCrVPw/Inl3CwKcrPokJNJfk4medmZ5OdkkZ+TSYfsTDrsf87ObNM/Mwg1vMxsLPAbEneVn+zuP22xPhf4A3AasBX4mruvCdbdBlwHNAE3uftLYdYahYrgeld5mcJLpL0wMzrnZdM5L5sBxQWHvV9tQ1MQakG47W1gT10jNcFjd20jNXUN1NT+dXnjrlpqqhupqW1kd13jIQOwpdysDPJzEuGWl51Bfk7WJwIuNzuDvKzgOTuT3Ky/PpcW5TP2xJ6f9T/PYQstvMwsE7gbOB+oAuab2TR3X5a02XXAdnc/xsyuAn4GfM3MhpCYHmUo0Bv4s5kd6+6fvGoacwvWbqe4Yy79uuZHXYqItHF52Zn0LMykZ2He536PusYm9tQ1saeukX0NTeytb2JvfSO1wet99U0ft+9/va8+sZzYppG99U3s2FvPhoYm6hqbqW3xvP/nCGcO6hbP8AJGAJXuvhrAzKYA44Hk8BoP3BG8ngrcZYmhOeOBKe5eB3wQzPc1gsSklaH48XNLWfbRrrDevlXLNuxi9KBijUYSkaMiNyuT3KxMuhbkhPL+7k59U3MixD7bSd5nFuYtHfoA65KWq4K2Vrdx90ZgJ9DtMPfFzK43swozq6iurk5h6eFranZyszI4Y6B+3yUi6cHMyM3KpHNedugjqGM9YMPdJwGTIDEZ5ZG81+2XDk1JTSIiEr4wz7zWA8lzfJQGba1uY2ZZQCGJgRuHs6+IiLRTYYbXfGCwmQ0wsxwSAzCmtdhmGjAheH0l8Kq7e9B+lZnlmtkAYDDwVoi1iohIjITWbejujWY2EXiJxFD5+919qZndCVS4+zTgPuChYEDGNhIBR7DdEyQGdzQCN6bbSEMREfn8zP2ILhW1GWZWDaw9wrcpBrakoJy2TMcYf+l+fKBjTBdHeoz93b2ktRVpE16pYGYV7l4edR1h0jHGX7ofH+gY00WYx6jZD0VEJHYUXiIiEjsKr0+aFHUBR4GOMf7S/fhAx5guQjtGXfMSEZHY0ZmXiIjEjsJLRERiR+EVMLOxZrbczCrN7Nao60k1M7vfzDab2btR1xIGM+trZjPNbJmZLTWzm6OuKdXMLM/M3jKzd4Jj/HHUNYXBzDLN7G0z+1PUtYTBzNaY2RIzW2RmFVHXEwYz62JmU83sfTN7z8xGpfwzdM3r47nHVpA09xhwdYu5x2LNzM4GaoA/uPuJUdeTambWC+jl7gvNrBOwALgszf4MDShw9xozywbeBG5297kRl5ZSZnYLUA50dvdLoq4n1cxsDVDu7mn7A2UzexB4w90nB7cHzHf3Han8DJ15JXw895i71wP75x5LG+7+OolbcKUld9/g7guD17uB92hlGp0484SaYDE7eKTVt08zKwW+DEyOuhb5fMysEDibxO3/cPf6VAcXKLz2O6z5wyQezKwMOAWYF3EpKRd0qS0CNgMz3D3djvHXwPeBkKcyjJQDL5vZAjO7PupiQjAAqAYeCLp/J5tZQao/ROElacXMOgJPAd9196M7NfZR4O5N7j6cxDRBI8wsbbqAzewSYLO7L4i6lpCNcfdTgYuAG4Mu/XSSBZwK3OPupwB7gJSPI1B4JWj+sDQQXAd6CnjE3Z+Oup4wBd0wM4GxEZeSSqOBccE1oSnAuWb2cLQlpZ67rw+eNwPPkLhskU6qgKqkXoGpJMIspRReCYcz95i0YcFghvuA99z9l1HXEwYzKzGzLsHrDiQGGL0faVEp5O63uXupu5eR+Df4qrtfG3FZKWVmBcGAIoKutAuAtBoB7O4bgXVmdlzQ9CUS01ulVGjzecXJgeYei7islDKzx4BzgGIzqwJud/f7oq0qpUYDXweWBNeEAH7o7tOjKynlegEPBqNjM4An3D0th5OnsR7AM4nvWmQBj7r7i9GWFIrvAI8EJwOrgW+m+gM0VF5ERGJH3YYiIhI7Ci8REYkdhZeIiMSOwktERGJH4SUiIrGj8BKJmJl1C+4wvsjMNprZ+uB1jZn9Pur6RNoiDZUXaUPM7A6gxt3/K+paRNoynXmJtFFmds7+Oa3M7A4ze9DM3jCztWZ2uZn9PJgX6sXg1liY2Wlm9pfgpq8vBVPFiKQdhZdIfAwCzgXGAQ8DM919GLAP+HIQYL8DrnT304D7gf+IqliRMOn2UCLx8YK7N5jZEhK3Mdt/W6ElQBlwHHAiMCO4/VAmsCGCOkVCp/ASiY86AHdvNrMG/+sF62YS/5YNWOruKZ9yXaStUbehSPpYDpSY2ShITBFjZkMjrkkkFAovkTTh7vXAlcDPzOwdYBFwZqRFiYREQ+VFRCR2dOYlIiKxo/ASEZHYUXiJiEjsKLxERCR2FF4iIhI7Ci8REYkdhZeIiMTO/w+qux/gxcYkEQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "B_y = 0\n", "Beta_y = 1\n", "B_z = 0\n", "Beta_z = 1\n", "K_xy = 0.1\n", "K_xz = 0.1\n", "K_yz = 0.3 # 1, 0.3, 0.1\n", "H = 2\n", "alpha_y = 1\n", "alpha_z = 1\n", "\n", "time_x_on = 0.5\n", "time_x_off = 2.5\n", "params_off = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,0,1);\n", "params_on = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,1,1);\n", "\n", "mf_off = lambda t,y : incoherent_type_4_and_model(t,y, params_off)\n", "mf_on = lambda t,y : incoherent_type_4_and_model(t,y, params_on)\n", "\n", "sol1 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(0,time_x_on), y0=(0,0), \n", " t_eval=np.arange(0,time_x_on,.01) )\n", "sol2 = scipy.integrate.solve_ivp( fun=mf_on, t_span=(time_x_on,time_x_off), y0=sol1.y[:,-1], \n", " t_eval=np.arange(time_x_on,time_x_off,.01) )\n", "sol3 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(time_x_off,6), y0=sol2.y[:,-1], \n", " t_eval=np.arange(time_x_off,6,.01) )\n", "plt.subplot( 2, 1, 1)\n", "plt.plot( [sol1.t[0],time_x_on,time_x_on,time_x_off,time_x_off,sol3.t[-1]], [0,0,1,1,0,0] )\n", "plt.ylabel( \"X\")\n", "plt.subplot( 2, 1, 2 )\n", "plt.plot(sol1.t, sol1.y[1,:].T); #Simulation before signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol2.t, sol2.y[1,:].T); #Simulation after signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol3.t, sol3.y[1,:].T); #Simulation after signal on\n", "plt.xlabel( 'Time')\n", "plt.ylabel( \"Z\")\n", "plt.tight_layout()\n", " " ] }, { "cell_type": "code", "execution_count": 41, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa8AAAEYCAYAAADrpHnMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqKElEQVR4nO3deXhU9d3+8fcnO0mABBLWAGFzAZQtgkrdq0WrwqO2LtWqj1bbn1pb+7TV9rm6WLva2k3bSgW1akXF2tKKWxV3QRIB2TEikIQlCUtC2LJ9fn/MwBNpQJaZnJzJ/bquuTJzziz3EeHOOec752vujoiISJgkBR1ARETkUKm8REQkdFReIiISOiovEREJHZWXiIiETkrQAWIlLy/PCwsLg44hIiIxUlJSUu3u+a2tS5jyKiwspLi4OOgYIiISI2a2Zn/r4nbY0MymmVmlmS3ez3ozs9+ZWamZvW9mY1qsu9rMPojero5XRhERCad4nvN6CJh4gPXnAkOjtxuAPwKYWTfg+8B4YBzwfTPLjWNOEREJmbgdNnT3182s8ABPmQT8xSOX+JhjZjlm1hs4HXjJ3TcDmNlLRErw8XhlDcpf567lHwsqgo4hR2jSqL5cMb5/0DFEOpQgRxv2BcpaPC6PLtvf8v9gZjeYWbGZFVdVVcUtaLz8Y0EFS9fXBh1DjsDS9bX6BUQkAKEesOHuU4ApAEVFRaG8SOOw3l144saTgo4hh+nS+98JOoJIhxTknlcF0K/F44Losv0tFxERAYItr5nAF6OjDk8Eatx9PfACcI6Z5UYHapwTXSYiIgLE8bChmT1OZPBFnpmVExlBmArg7n8CZgHnAaXADuDa6LrNZvYjYF70re7cM3hDREQE4jva8PJPWO/ATftZNw2YFo9cIiISfrq2oYiIhI7KS0REQkflJSIioaPyEhGR0FF5iYhI6Ki8REQkdFReIiISOiovEREJHZWXiIiEjspLRERCR+UlIiKho/ISEZHQUXmJiEjoqLxERCR0VF4iIhI6cS0vM5toZivMrNTMbm9l/a/NbEH0ttLMtrZY19Ri3cx45hQRkXCJ50zKycB9wNlAOTDPzGa6+9I9z3H3r7d4/i3A6BZvsdPdR8Urn4iIhFc897zGAaXuvsrd64HpwKQDPP9y4PE45hERkQQRz/LqC5S1eFweXfYfzGwAMBB4pcXiDDMrNrM5ZjZ5P6+7Ifqc4qqqqhjFFhGR9q69DNi4DJjh7k0tlg1w9yLgCuA3ZjZ43xe5+xR3L3L3ovz8/LbKKiIiAYtneVUA/Vo8Logua81l7HPI0N0roj9XAa/y8fNhIiLSgcWzvOYBQ81soJmlESmo/xg1aGbHALnAOy2W5ZpZevR+HjABWLrva0VEpGOK22hDd280s5uBF4BkYJq7LzGzO4Fid99TZJcB093dW7z8WOB+M2smUrA/azlKUUREOra4lReAu88CZu2z7Hv7PP5BK697GzguntlERCS82suADRERkYOm8hIRkdBReYmISOiovEREJHRUXiIiEjoqLxERCR2Vl4iIhI7KS0REQkflJSIioaPyEhGR0FF5iYhI6Ki8REQkdFReIiISOiovEREJnbiWl5lNNLMVZlZqZre3sv4aM6syswXR2/Ut1l1tZh9Eb1fHM6eIiIRL3ObzMrNk4D7gbKAcmGdmM1uZVPIJd795n9d2A74PFAEOlERfuyVeeUVEJDz2u+dlZv0PsO6Ug3jvcUCpu69y93pgOjDpIHN9BnjJ3TdHC+slYOJBvlZERBLcgQ4bvmpm34ruQQFgZj3N7FHg1wfx3n2BshaPy6PL9nWxmb1vZjPMrN+hvNbMbjCzYjMrrqqqOohIIiKSCA5UXmOBwcACMzvTzG4F3gXeIbJXFQv/BArd/Xgie1cPH8qL3X2Kuxe5e1F+fn6MIomISHu333Ne0cN1N0ZL69/AOuBEdy8/yPeuAPq1eFwQXdbyMza1ePgA8IsWrz19n9e+epCfKyIiCe5A57xyzOx+4Foi55tmAM+Z2ZkH+d7zgKFmNtDM0oDLgJn7fEbvFg8vBJZF778AnGNmuWaWC5wTXSYiInLA0YbvAX8AbnL3RuBFMxsF/MHM1rj75Qd6Y3dvNLObiZROMjDN3ZeY2Z1AsbvPBL5qZhcCjcBm4Jroazeb2Y+IFCDAne6++bC3UkREEsqByuvUfQ8RuvsC4GQz+9LBvLm7zwJm7bPsey3u3wHcsZ/XTgOmHczniIhIx7Lfw4YHOrfl7n+OTxwREZFPpstDiYhI6Ki8REQkdFReIiISOiovEREJHZWXiIiEjspLRERCR+UlIiKho/ISEZHQUXmJiEjoqLxERCR0VF4iIhI6Ki8REQkdlZeIiISOyktEREInruVlZhPNbIWZlZrZ7a2sv83MlprZ+2b2spkNaLGuycwWRG8z932tiIh0XAeajPKImFkycB9wNlAOzDOzme6+tMXT5gNF7r7DzL4C/AK4NLpup7uPilc+EREJr3jueY0DSt19lbvXA9OBSS2f4O6z3X1H9OEcoCCOeUREJEHEs7z6AmUtHpdHl+3PdcBzLR5nmFmxmc0xs8mtvcDMbog+p7iqquqIA4uISDjE7bDhoTCzK4Ei4LQWiwe4e4WZDQJeMbNF7v5hy9e5+xRgCkBRUZG3WWAREQlUPPe8KoB+LR4XRJd9jJl9GvgucKG7796z3N0roj9XAa8Co+OYVUREQiSe5TUPGGpmA80sDbgM+NioQTMbDdxPpLgqWyzPNbP06P08YALQcqCHiIh0YHE7bOjujWZ2M/ACkAxMc/clZnYnUOzuM4G7gWzgKTMDWOvuFwLHAvebWTORgv3ZPqMURUSkA4vrOS93nwXM2mfZ91rc//R+Xvc2cFw8s4mISHjpChsiIhI6Ki8REQkdlZeIiISOyktEREJH5SUiIqGj8hIRkdBReYmISOiovEREJHRUXiIiEjoqLxERCR2Vl4iIhI7KS0REQkflJSIioaPyEhGR0IlreZnZRDNbYWalZnZ7K+vTzeyJ6Pq5ZlbYYt0d0eUrzOwz8cwpIiLhErfyMrNk4D7gXGAYcLmZDdvnadcBW9x9CPBr4OfR1w4jMvPycGAi8Ifo+4mIiMR1MspxQKm7rwIws+nAJKDljMiTgB9E788A7rXIlMqTgOnuvhv4yMxKo+/3TrzC/vCfS1i6rjZeb9+qpetrGda7S5t+psTe0vW1XHp/3P7XFAmlYX268P0Lhsft/eNZXn2BshaPy4Hx+3uOuzeaWQ3QPbp8zj6v7bvvB5jZDcANAP37949Z8LYyrHcXJo36j82SENGfn0gw4llecefuU4ApAEVFRX4k7xXP3xAkcV0xvj9XjA/fL04iYRfPARsVQL8Wjwuiy1p9jpmlAF2BTQf5WhER6aDiWV7zgKFmNtDM0ogMwJi5z3NmAldH718CvOLuHl1+WXQ04kBgKPBuHLOKiEiIxO2wYfQc1s3AC0AyMM3dl5jZnUCxu88EpgKPRAdkbCZScESf9ySRwR2NwE3u3hSvrCIiEi4W2dEJPzOrAtYc4dvkAdUxiNOeaRvDL9G3D7SNieJIt3GAu+e3tiJhyisWzKzY3YuCzhFP2sbwS/TtA21joojnNuryUCIiEjoqLxERCR2V18dNCTpAG9A2hl+ibx9oGxNF3LZR57xERCR0tOclIiKho/ISEZHQUXlFfdLcY2FnZtPMrNLMFgedJR7MrJ+ZzTazpWa2xMxuDTpTrJlZhpm9a2YLo9v4w6AzxYOZJZvZfDP7V9BZ4sHMVpvZIjNbYGbFQeeJBzPLMbMZZrbczJaZ2Ukx/wyd89o799hK4GwiV7CfB1zu7ksP+MIQMbNTgTrgL+4+Iug8sWZmvYHe7v6emXUGSoDJCfZnaECWu9eZWSrwJnCru8/5hJeGipndBhQBXdz9/KDzxJqZrQaK3D1hv6BsZg8Db7j7A9HLA2a6+9ZYfob2vCL2zj3m7vXAnrnHEoa7v07kElwJyd3Xu/t70fvbgGW0Mo1OmHlEXfRhavSWUL99mlkB8FnggaCzyOExs67AqUQu/4e718e6uEDltUdrc48l1D98HYmZFQKjgbkBR4m56CG1BUAl8JK7J9o2/gb4FtAccI54cuBFMyuJzkmYaAYCVcCD0cO/D5hZVqw/ROUlCcXMsoGnga+5e9tOjd0G3L3J3UcRmSZonJklzCFgMzsfqHT3kqCzxNmn3H0McC5wU/SQfiJJAcYAf3T30cB2IObjCFReEZo/LAFEzwM9DTzm7n8LOk88RQ/DzAYmBhwlliYAF0bPCU0HzjSzR4ONFHvuXhH9WQk8Q+S0RSIpB8pbHBWYQaTMYkrlFXEwc49JOxYdzDAVWObu9wSdJx7MLN/McqL3OxEZYLQ80FAx5O53uHuBuxcS+Tv4irtfGXCsmDKzrOiAIqKH0s4BEmoEsLtvAMrM7OjoorOITG8VU3GbzytM9jf3WMCxYsrMHgdOB/LMrBz4vrtPDTZVTE0ArgIWRc8JAXzH3WcFFynmegMPR0fHJgFPuntCDidPYD2BZyK/a5EC/NXdnw82UlzcAjwW3RlYBVwb6w/QUHkREQkdHTYUEZHQUXmJiEjoqLxERCR0VF4iIhI6Ki8REQkdlZeIiISOyktEREJH5SUiIqGj8hIRkdBReYmISOiovEREJHQS5sK8eXl5XlhYGHQMERGJkZKSkmp3z29tXcKUV2FhIcXFxUHHEBGRGDGzNftbp8OGIiISOgmz5yWJy91pbjFzT5JBdD4kEemgVF7SLuysb2JRRQ2LKmpYXb2d1Zu2s75mF1t3NFCzs56Gpv9rr7SUJHIzU8nNTCMvO53CvEwG5mUzKD+LEX26kt85PcAtEZG2oPKSQDQ3OwvKt/LKskpeW1nF0vW1NEV3rzpnpDAwL4uhPbLJyUwjJzOVjJRkzMAddtQ3smVHPVt2NFBZu4t/LFjHtl2Ne9+7X7dOjO6XywkDu3H6Ufn065YZ1GaKSJyovKRNlVZu46nicv42v4KqbbtJTjLG9s/lK6cNZnT/HI4vyCEvO+2QDgu6O5u311NaWcf75TXML9vCux9tZubCdQAMys/ijKN7cN5xvRnTP0eHHEUSgLn7Jz8rBIqKilyjDdun5mbn5eWVTHn9Q+at3kJKknHGMT04//jenHZUPjmZaTH/THdnVfV2Xl1RxasrKpm7ajP1Tc0U5HbigpF9uHhMX4b06BzzzxWR2DGzEncvanWdykvipbnZ+ef767j3lVI+qKyjILcTXzxpAP81uqDNz0vV7mrgxSUb+efCdbxZWk1Ts3PioG58YfwAPjO8F2kpGngr0t6ovKTNvfPhJn48aymLK2o5pldnvnL6YD57XG9SkoMvieq63TxVXM5f311D2ead5GWnccX4AVxzciHdsmK/Fygih0flJW2mctsufjBzCbMWbaBP1wy+OfFoJo3sS1JS+zvP1NzsvFFazSPvrObfyyrJSE3i0qJ+XH/KIA3yEGkHVF4Sd+7OUyXl3PWvpexqbOaWM4bwpVMHkZGaHHS0g1JauY37X1vF3xdU0OwweVRfbj1rKP27q8REgqLykrjasr2eb85YyL+XVTKusBs/vfg4BudnBx3rsKyv2ckDb3zEo3PW0NTsXHpCP245cyi9umYEHU2kw1F5Sdy8+9Fmbp0+n0119dx+7jFcc3JhuzxEeKg21Ozi3tkf8MS8MsyMq04cwC1nDonLyEgRaZ3KS2LO3XnwrdXc9exS+nfL5N4rxjCib9egY8Vc2eYd/PblD/jbe+V0zkjl1rOGctVJA0htBwNPRBKdyktiqr6xme/9YzHT55VxzrCe3HPpKLLTE/v77ss31HLXv5bxZmk1g/Ky+M55x3LWsT30hWeRODpQeenXRzkkNTsbuGrqXKbPK+PmM4bwpyvHJnxxARzTqwuPXDeOadcUgcH1fynmyqlzWbGhls3b64OOJ9LhaM9LDlrltl18ceq7fFhVx92XjGTy6L5BRwpEQ1Mzj81Zw6///QG1Oxvo0SWd2f9zOplpiV/iIm1Je15yxMo27+Bzf3qHtZt3MO2aEzpscQGkJidxzYSBzP6f0zlneE8qa3fzo38tDTqWSIei8pJPtLp6O5f86W227mjg0evHc8rQVmfl7nC6ZaVx/1VFfPn0wTz+bhnPLVofdCSRDkPlJQdUtnkHV/x5DvWNzTxx44mM6Z8bdKR257azj2JkQVdu/9si1m3dGXQckQ5B5SX7tW7rTi7/8xy21zfx6PXjOaZXl6AjtUupyUn87vLRNDY187UnFuydl0xE4kflJa3avL2eKx+YS82OBh65bhzD+yTed7hiaUD3LH40eQTvfrSZ+2aXBh1HJOGpvOQ/7Kxv4vqH51G+dSdTrzmB4wtygo4UCheNKWDyqD785t8rebu0Oug4IglN5SUf09Ts3Dp9PvPLtvLbS0cxbmC3oCOFyo//K3Jdx1sen8/6Gp3/EokXlZd8zF3PLuXFpRv5wQXDOfe43kHHCZ2s9BT+eOVYdjU0cdNj71Hf2Bx0JJGEFEh5mdlEM1thZqVmdvsBnnexmbmZtfolNYmtJ4vLePCt1fz3hIFcfXJh0HFCa0iPbH5xyUjeW7uVn8xaFnQckYTU5uVlZsnAfcC5wDDgcjMb1srzOgO3AnPbNmHHNH/tFv73mcVMGNKd75x3TNBxQu+zx/fmuk8N5KG3VzNz4bqg44gknCD2vMYBpe6+yt3rgenApFae9yPg58CutgzXEVXW7uLLj5bQs2s6914+hhRdMT0mbj/3GE4ozOXbM95ncUVN0HFEEkoQ/0r1BcpaPC6PLtvLzMYA/dz92QO9kZndYGbFZlZcVVUV+6QdQENTM19+tITanY1MuaqI3CzNVxUrqclJ/OELY8nNTOVLfymmcpt+DxOJlXb3K7aZJQH3AN/4pOe6+xR3L3L3ovx8XbLocPzqxZW8t3Yrv7jkeI7trS8hx1p+53T+fHURW3c0cMNfStjV0BR0JJGEEER5VQD9WjwuiC7bozMwAnjVzFYDJwIzNWgj9l5fWcWfXvuQy8f154KRfYKOk7CG9+nKry8dyYKyrdzxt0UkykwOIkEKorzmAUPNbKCZpQGXATP3rHT3GnfPc/dCdy8E5gAXurvmO4mhqm27ue3JhRzVM5vvnf8f42UkxiaO6M3/nHMUz8yv4A+vfhh0HJHQa/MJiNy90cxuBl4AkoFp7r7EzO4Eit195oHfQY5Uc7Nz25ML2LargceuH0+ntOSgI3UIN50xhA8q67j7hRX0zenUoaeVETlSgcye5+6zgFn7LPvefp57eltk6kgefHs1b3xQzV2TR3B0r85Bx+kwzIxfXHI8G2t38c0ZC8nvnM6EIXlBxxIJpXY3YEPi66Pq7dz9wnLOPKYHXxjfP+g4HU56SjL3X1XE4PxsbnykhKXraoOOJBJKKq8OpLnZ+daMhaQlJ/HTi47DzIKO1CF17ZTKQ9eOo0tGCtc8+C7lW3YEHUkkdFReHciDb69m3uotfO+C4fTskhF0nA6tV9cMHvrvcexqaOKLU9+latvuoCOJhIrKq4Noebjw4jEaKNAeHNWzM9OuOYH1Nbu4aupctu6oDzqSSGiovDoAd+fbT7+vw4XtUFFhNx64uohV1du5etq7bNvVEHQkkVBQeXUAT79XwbsfbeY75x2rw4Xt0IQhefzxC2NYsq6W/35oHjvqG4OOJNLuqbwS3NYd9fx01jLG9M/h80X9PvkFEoizju3Jby4bRcmaLVz/cLEKTOQTqLwS3N0vrGDLjnrumnwcSUk6XNienX98H375uZHMWbWJax6cR91uFZjI/qi8EtiCsq389d21XHPyQIb10UV3w+CiMQX89rLRlKzZwhenzqVW58BEWqXySlBNzc53n1lEfnY6Xz97aNBx5BBcMLIP910xhkUVNVz5gEYhirRG5ZWgniwuY8m6Wv73/GF0zkgNOo4cookjevGnK8eyfP02Lv/zXKrr9D0wkZZUXgmobncjv3pxBUUDcrng+N5Bx5HDdNaxPXng6iI+qq7j4j++zZpN24OOJNJuqLwS0J9e/ZDqunq++9lj9Z2ukDv1qHz++qUTqd3ZwMV/fJtF5TVBRxJpF1ReCaZi607+/MYqJo3qw+j+uUHHkRgY0z+XGV85mfSUZC6d8g6vrawKOpJI4FReCebu55cD8K2JxwScRGJpcH42z/y/kxnQPYvrHprH0yXlQUcSCdRhlZeZfeYA6z53+HHkSCws28rfF6zj+lMG0jenU9BxJMZ6dMngyRtPZPygbnzjqYX8/PnlNDd70LFEAnG4e16zzGy2mbV2hdc7jiSQHL6fPbecvOw0vnL6kKCjSJx0zohMp3LF+P788dUPueGRYn2ZWTqkwy2v94G/AnPM7JJ91mmEQADeLq3mnVWb+H+nDyE7PZAJsqWNpCYn8ePJI7hz0nBmr6ji4j+8TdlmzQkmHcvhlpe7+5+Bs4Bvm9mDZpa5Z11sosnBcnfufnEFvbtmcIVmR+4QzIwvnlTIw9eOY33NTi68903eKq0OOpZImzmiARvuvhI4CdgIzDez8TFJJYdk9opK5q/dylfPGkpGanLQcaQNfWpoHv+4+VN0z07nqqlz+f3LH+g8mHQIh1teew8Nunuju98O3Ag8DuhaRG2oudn55QsrGdA9k0vGFgQdRwIwMC+Lf9w0gQtG9uFXL63k2ofmsXm7Liklie1wy+uH+y5w91eBscCPjySQHJrnl2xg6fpavvbpoaQm65sPHVVWegq/uXQUd00ewTsfbuL8373Be2u3BB1LJG4O6187d//7fpZvcfefHVEiOWhNzc49L61kaI9sLhzZ2sBP6UjMjCtPHMCMr5xEUpLxuT+9w+9f/oAmHUaUBKRf1UPs+cUbKK2s49ZPDyVZc3VJ1PEFOTx7yymcd1xvfvXSSi69/x2NRpSEo/IKKXfnvtmlDMrP4twRuviufFzXzFR+d9kofnPpKFZs2Ma5v32Dp0vKcddemCQGlVdIvbqiiqXra/nKaYO11yWtMjMmj+7LrFtPYVjvLnzjqYXc+EgJG2t3BR1N5IgFUl5mNtHMVphZqZnd3sr628xsqZm9b2Yvm9mAIHK2V+7OvbNL6ZvTicmjda5LDqxft0wev+FE7jj3GF5bWcWn73mNJ+at1V6YhFqbl5eZJQP3AecCw4DLzWzYPk+bDxS5+/HADOAXbZuyfZuzajMla7bw5dMGaYShHJTkJOPG0wbz/NdO5djeXfj204u4cupc1m7SuTAJpyD+5RsHlLr7KnevB6YDk1o+wd1nu/uev1VzAH2BqYX7ZpeSl53O54r6BR1FQmZgXhbTv3Qid00ewcKyGj7zm9e5/7UPqW9sDjqayCEJorz6AmUtHpdHl+3PdcBzra0wsxvMrNjMiquqOsYcRwvLtvJmaTVfOmWgrqYhhyUpKTKk/sWvn8qEIXn89LnlTPzt67yuecIkRNr1MSczuxIoAu5ubb27T3H3Incvys/Pb9twAZnyxio6Z6TwhRN1GlCOTJ+cTjxwdREPXnMCzc3OF6e9y42PFGtYvYRCEOVVAbQ83lUQXfYxZvZp4LvAhe6+u42ytWtlm3fw3KL1XDG+v64cLzFzxjE9eOHrp/LNzxzN6yur+fQ9r3HPiys01Yq0a0GU1zxgqJkNNLM04DJgZssnmNlo4H4ixVUZQMZ26eG3V5NkxjUnFwYdRRJMekoyN50xhJe/cRrnDO/F714p5bRfzObht1frfJi0S21eXu7eCNwMvAAsA5509yVmdqeZXRh92t1ANvCUmS0ws5n7ebsOY9uuBqbPK+O843rTu6tmSZb46JPTid9fPpp/3DSBo3p25vszl3D2r1/jnwvX6Wr10q5YonzXo6ioyIuLi4OOETcPvLGKu55dxsybJ3B8QU7QcaQDcHdeXVnFz59bzvIN2xjepwtfPWso5wzriZm+GC/xZ2Yl7l7U2rp2PWBDIhqbmnnwrdWMK+ym4pI2Y2accXQPnv3qKfzqcyPZvruRGx8p4bzfvclzi9ZrT0wCpfIKgReWbKRi606uO2Vg0FGkA0pOMi4eW8C/bzuNez4/kt0NTXzlsfc473dv8M+F62hs0jkxaXsashYCU99cxYDumXz62J5BR5EOLCU5iYvGFDBpVF/+9f46fvfyB9zy+Hz65nTi2gmFXDZOo2Cl7WjPq51bXFHDe2u38sWTCnUBXmkXkpOMSaP68tLXT2PKVWPpm9OJu55dxkk/eZmfzFrGuq07g44oHYB+TWrnHpu7hozUJC4ZoytkSfuSlGScM7wX5wzvxcKyrfz5jVVMffMjpr35EZ8Z3osvjO/PSYO7a3CHxIXKqx2r2dnA3+evY9LIvnTNTA06jsh+jeyXw71XjKF8yw4eems1T5WU8+yi9QzKy+KK8f25ZGwBOZlpQceUBKLDhu3Y0yXl7Gxo4qqTdCkoCYeC3Ez+9/xhzP3OWdzz+ZHkZqVx17PLGPeTl7ntiQW8XVqtUYoSE9rzaqfcnUfnrGFUvxxG9O0adByRQ5KRmsxFYwq4aEwBy9bX8te5a/n7/Ar+Nr+C3l0zmDy6LxeN7svQnp2DjiohpS8pt1NvlVbzhQfmcs/nR3KRzndJAtjV0MRLSzfyzPwKXltZRVOzc1zfrvzX6L6cf3xvenTJCDqitDMH+pKyyqud+vIjJcz9aBPv3HGWpj6RhFO1bTczF67jmfnlLK6oxQxOGNCNc4/rxcQRvXQJNAFUXqGzoWYXE37+CtefMpA7zj026DgicVVauY1ZizYwa9F6lm/YBsDo/jmcN6I3nxnei/7dMwNOKEFReYXMva98wC9fXMlr3zydAd2zgo4j0mZWVdXx3OINPLd4PYsragEYlJ/FmUf34MxjelBU2I20FI0z6yhUXiHS3Oyc/stX6ZOTwfQbTgo6jkhg1m7awcvLN/LK8krmrtpMfVMz2ekpnDI0jzOO6cGpQ/Pp1VXnyRLZgcpLow3bmTkfbWLt5h18/eyhQUcRCVT/7plcO2Eg104YyPbdjbxVWs3sFZW8sryS5xZvACJ7ZScP7s6EwXmcNLi7vkvWgai82pkn55XROSOFc0f0DjqKSLuRlZ6y92oe7s6y9dt4+8Nq3iqt5pn3Knh0zlrMYFjvLkwYkkfRgFzGDsile3Z60NElTlRe7UjNzgaeW7yBzxUVaIShyH6YGcP6dGFYny5cf8ogGpqaeb98K2+VbuKt0moeems1U15fBcCgvCzGDMilaEAuRYW5DMrLJknXCE0IKq92ZOaCCnY3NnNpUf+go4iERmpyEmMHdGPsgG589ayh7GpoYlFFDSVrtlC8egsvL9vIjJJyAHIyUxlZkMNxfbsyom9XjivoSp+uGbr+YgipvNqRJ4rLOLZ3F0b07RJ0FJHQykhN5oTCbpxQ2A1Oi1ytZlX1dkrWbKFk9RYWlm/lzdJqmqKXqeqelRYpsmihDevdhYLcTtpDa+dUXu3EknU1LK6o5YcXDtdvgSIxZGYMzs9mcH42ny/qB0Su9rFsfS2LKmpYVF7DooqajxVaZloyQ3t25uie2RzVszPH9OrCUb2yyc9O19/PdkLl1U48VVxOWkoSk0b1CTqKSMLLSE1mdP9cRvfP3btsV0MTS9fXsmLDNlZs2MbKjdt4eVklTxaX731ObmYqQ3t2ZlBeFgOjt0H5WfTrlkl6is5TtyWVVzvQ0NTMPxeu4+xje2qor0hAMlKTGdM/lzEtCg2gum43KzdsY8XGSKGt3FjHS0s3sml7/d7nJFnkivp7Cq2weyb9umVSkJtJQW4nsjTDdMzpv2g78OYH1WzaXs/k0X2DjiIi+8jLTidvSDonD8n72PKaHQ18tGk7H1XX8VH1Dj6qjtwvXr2Z7fVNH3tubmbq3iKL3DL3/uzVNYMuGSk6HHmIVF7twN/mV5CbmcppR+UHHUVEDlLXzFRGZeYwql/Ox5a7O9V19ZRv2UH5lp3RW+T+io3beHl5JfWNzR97TafUZHp2Sadnlwx6dc2gZ5fIrVeXDHp1jSzv0TlDl8ZqQeUVsG27GnhxyQY+X9RP/2OKJAAzI79zOvmd0z92Tm2P5manevvuvcW2sWYXG2p3sTF6e2/tFjbW7v6PggPokpFCXnY63bPT6J6152ca3Vssy8uOPM7plJrQIyZVXgF7fvEGdjc265ChSAeRlGT06BzZk9r3/Noe7s7WHQ1sqI0UW2XtLjbW7mbz9nqq63azqa6eVdV1zFtdz+Yd9bR2iVoz6JKRStdOqeRkRn7uuX38cdreZV06pZKdnkJ2egrJ7bz4VF4B+/uCCgZ0z2RM/5ygo4hIO2Fm5GalkZuVxrG9D/y9z6ZmZ8uOejbV1bOpbjfV2+vZXLebTdvrqdnZQM3OBrbuiPws37Jz77I9XwvYn8y0ZLLSU+icnkJ2RgpZaZGfndNTyIou21N02dFlWenJZKalkJmWTE5malznZQukvMxsIvBbIBl4wN1/ts/6dOAvwFhgE3Cpu69u65zxtr5mJ29/uImvnjlUJ2tF5LAkJ1lkUEl2OtD5oF7j7tTtbtxbbLU7G9gaLbXtuxup291I3a7oz+ht++5GyjbvYHv9/61raNp/AY4f2I0nbozfzBhtXl5mlgzcB5wNlAPzzGymuy9t8bTrgC3uPsTMLgN+Dlza1lnjbeaCdbijQ4Yi0qbMjM4ZqXTOSKWg9SOXB2V3YxN1uxrZvruJbbsb2FHfFLntbqRLp9TYBW5FEHte44BSd18FYGbTgUlAy/KaBPwgen8GcK+Zmcdx8rEf/nMJS9fVxuvtW7Vy4zZG989hYJ4mnBSR8ElPSSY9O5nu2W3/2UEMb+sLlLV4XB5d1upz3L0RqAG67/tGZnaDmRWbWXFVVVWc4sZHY3MzKUnGRdrrEhE5ZKEesOHuU4ApEJlJ+Uje6/sXDI9JpkOVKDNZi4i0pSD2vCqAfi0eF0SXtfocM0sBuhIZuJFwNFBDROTQBVFe84ChZjbQzNKAy4CZ+zxnJnB19P4lwCvxPN8lIiLh0uaHDd290cxuBl4gMlR+mrsvMbM7gWJ3nwlMBR4xs1JgM5GCExERAcASZYfGzKqANUf4NnlAdQzitGfaxvBL9O0DbWOiONJtHODurV70NWHKKxbMrNjdi4LOEU/axvBL9O0DbWOiiOc26kqwIiISOiovEREJHZXXx00JOkAb0DaGX6JvH2gbE0XctlHnvEREJHS05yUiIqGj8hIRkdBReUWZ2UQzW2FmpWZ2e9B5Ys3MpplZpZktDjpLPJhZPzObbWZLzWyJmd0adKZYM7MMM3vXzBZGt/GHQWeKBzNLNrP5ZvavoLPEg5mtNrNFZrbAzIqDzhMPZpZjZjPMbLmZLTOzmE/spXNe7J1jbCUt5hgDLt9njrFQM7NTgTrgL+4+Iug8sWZmvYHe7v6emXUGSoDJCfZnaECWu9eZWSrwJnCru88JOFpMmdltQBHQxd3PDzpPrJnZaqDI3RP2C8pm9jDwhrs/EL0MYKa7b43lZ2jPK2LvHGPuXg/smWMsYbj760QutZWQ3H29u78Xvb8NWMZ/TrUTah5RF32YGr0l1G+fZlYAfBZ4IOgscnjMrCtwKpHL/OHu9bEuLlB57XEwc4xJSJhZITAamBtwlJiLHlJbAFQCL7l7om3jb4BvAc0B54gnB140sxIzuyHoMHEwEKgCHowe/n3AzGI+467KSxKKmWUDTwNfc/e2nRq7Dbh7k7uPIjKV0DgzS5hDwGZ2PlDp7iVBZ4mzT7n7GOBc4KboIf1EkgKMAf7o7qOB7UDMxxGovCIOZo4xaeei54GeBh5z978FnSeeoodhZgMTA44SSxOAC6PnhKYDZ5rZo8FGij13r4j+rASeIXLaIpGUA+UtjgrMIFJmMaXyijiYOcakHYsOZpgKLHP3e4LOEw9mlm9mOdH7nYgMMFoeaKgYcvc73L3A3QuJ/B18xd2vDDhWTJlZVnRAEdFDaecACTUC2N03AGVmdnR00VlAzAdOtfl8Xu3R/uYYCzhWTJnZ48DpQJ6ZlQPfd/epwaaKqQnAVcCi6DkhgO+4+6zgIsVcb+Dh6OjYJOBJd0/I4eQJrCfwTHQG9RTgr+7+fLCR4uIW4LHozsAq4NpYf4CGyouISOjosKGIiISOyktEREJH5SUiIqGj8hIRkdBReYmISOiovEQCZmbdo1cYX2BmG8ysInq/zsz+EHQ+kfZIQ+VF2hEz+wFQ5+6/DDqLSHumPS+RdsrMTt8zp5WZ/cDMHjazN8xsjZldZGa/iM4L9Xz00liY2Vgzey160dcXolPFiCQclZdIeAwGzgQuBB4FZrv7ccBO4LPRAvs9cIm7jwWmAT8OKqxIPOnyUCLh8Zy7N5jZIiKXMdtzWaFFQCFwNDACeCl6+aFkYH0AOUXiTuUlEh67Ady92cwa/P9OWDcT+btswBJ3j/mU6yLtjQ4biiSOFUC+mZ0EkSlizGx4wJlE4kLlJZIg3L0euAT4uZktBBYAJwcaSiRONFReRERCR3teIiISOiovEREJHZWXiIiEjspLRERCR+UlIiKho/ISEZHQUXmJiEjo/H+iQXVPAfNJ6gAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "B_y = 0\n", "Beta_y = 1\n", "B_z = 0\n", "Beta_z = 1\n", "K_xy = 0.1\n", "K_xz = 0.1\n", "K_yz = 0.1 # 1, 0.3, 0.1\n", "H = 2\n", "alpha_y = 1\n", "alpha_z = 1\n", "\n", "time_x_on = 0.5\n", "time_x_off = 2.5\n", "params_off = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,0,1);\n", "params_on = (B_y,Beta_y,B_z,Beta_z,K_xy,K_xz,K_yz,H,alpha_y,alpha_z,1,1);\n", "\n", "mf_off = lambda t,y : incoherent_type_4_and_model(t,y, params_off)\n", "mf_on = lambda t,y : incoherent_type_4_and_model(t,y, params_on)\n", "\n", "sol1 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(0,time_x_on), y0=(0,0), \n", " t_eval=np.arange(0,time_x_on,.01) )\n", "sol2 = scipy.integrate.solve_ivp( fun=mf_on, t_span=(time_x_on,time_x_off), y0=sol1.y[:,-1], \n", " t_eval=np.arange(time_x_on,time_x_off,.01) )\n", "sol3 = scipy.integrate.solve_ivp( fun=mf_off, t_span=(time_x_off,6), y0=sol2.y[:,-1], \n", " t_eval=np.arange(time_x_off,6,.01) )\n", "plt.subplot( 2, 1, 1)\n", "plt.plot( [sol1.t[0],time_x_on,time_x_on,time_x_off,time_x_off,sol3.t[-1]], [0,0,1,1,0,0] )\n", "plt.ylabel( \"X\")\n", "plt.subplot( 2, 1, 2 )\n", "plt.plot(sol1.t, sol1.y[1,:].T); #Simulation before signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol2.t, sol2.y[1,:].T); #Simulation after signal on\n", "plt.gca().set_prop_cycle(None) # this resets the color cycle\n", "plt.plot(sol3.t, sol3.y[1,:].T); #Simulation after signal on\n", "plt.xlabel( 'Time')\n", "plt.ylabel( \"Z\")\n", "plt.tight_layout()\n", " " ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }