chap6 updated.

Ziaeemehr · Aug 14, 2024 · 01b12cb · 01b12cb
1 parent 35e1da4
commit 01b12cb
Showing 1 changed file with 118 additions and 1 deletion.
diff --git a/docs/examples/chap_06.ipynb b/docs/examples/chap_06.ipynb
@@ -259,9 +259,44 @@
     "#### Stabiliby of steady states"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The process includes the calculation of the Jacobian matrix, finding the eigenvalues of the Jacobian matrix at the equilibrium point, and classifying the equilibrium point based on the nature of the eigenvalues.\n",
+    "\n",
+    "We'll use the `sympy` library for symbolic differentiation (to compute the Jacobian) and `numpy` for numerical linear algebra (to find eigenvalues).\n",
+    "\n",
+    "Here’s a step-by-step Python code that implements the process:\n",
+    "\n",
+    "### Explanation:\n",
+    "\n",
+    "1. **Step 1 (Equilibrium Points):** \n",
+    "   - We define a system of equations $ \\frac{dX}{dt} = F(X) $. \n",
+    "   - The equilibrium points are found by solving $ F(X) = 0 $.\n",
+    "\n",
+    "2. **Step 2 (Jacobian Matrix):**\n",
+    "   - The Jacobian matrix is calculated symbolically using `sympy.jacobian`.\n",
+    "\n",
+    "3. **Step 3 (Eigenvalue Analysis and Classification):**\n",
+    "   - The Jacobian matrix at each equilibrium point is substituted into and evaluated numerically.\n",
+    "   - The eigenvalues are calculated using `numpy.linalg.eigvals`.\n",
+    "   - The type of equilibrium point is classified based on the eigenvalues:\n",
+    "     - **Node (Stable/Unstable):** All real and either positive (unstable) or negative (stable).\n",
+    "     - **Saddle Point:** Mixed signs among real eigenvalues.\n",
+    "     - **Spiral (Stable/Unstable):** Complex eigenvalues with negative or positive real parts.\n",
+    "     - **Center:** Purely imaginary eigenvalues (neutral stability).\n",
+    "\n",
+    "### Customizing the Code:\n",
+    "\n",
+    "1. **System of Equations:** Modify `F` for your system.\n",
+    "2. **Variables:** If your system has more variables, extend the variables `x1`, `x2`, etc.\n",
+    "3. **Jacobian Calculation:** The Jacobian calculation automatically scales with the number of variables in the system."
+   ]
+  },
   {
    "cell_type": "code",
-   "execution_count": 54,
+   "execution_count": 2,
    "metadata": {},
    "outputs": [
     {
@@ -411,6 +446,88 @@
     "    print(f\"Equilibrium point {X_eq} is classified as: {classification}\")\n"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": 8,
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "'Equilibrium point:'"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAL0AAAAzCAYAAADCbQmYAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAJ00lEQVR4Ae2c7XXVOBCGbzgUwEIFGzqApYINHRCoYKEDOPwK/3Kgg2UrWKADoIOFDqCDhXSw+z6OlWvLkqWxZefmRjrH8ddoZvTOaDSSdXNwcnJyaxMor169Ogs8ro8qAlcGAflw0LdvqgVfdBx6LXmr+2fes+SthHwU0QudvyaJK8GmNcp7QXGs6xpkyvsE2B55bD8dKNJ/08NnAv2T99J0q/oI+KgzHeZaF2FwT0ev4+ueqHNb5+9dcKDVPdjd13V1/C44ha+F71OxPCbSzy5ihtG+63ztHb4F87OwwMmd43NNuX9+2v4V3VcdL/SEEffu9k29WgqBG3MZy2D0nsPWcHPZ7Uv9H2oIUZsojsN/0BGN5MKO9x90JnjUsjACsyK9jMRc4LWOQQRbWO9dZ0/0PrYoKXrmQt90PNVRR0wLeEbauZGeieupjNTLU406VPItAqQ5r4WnS4e2b+pVMQQmOz0RSVowMXtTTJtrzkhYkuYQQBg9a1kIgTnpDYapw3DEMG1QIGLf0UEayIjoJraRWs1joj2rYET8OoKOITXx3aRIL2M8lzwMejpR7r5XA5t3wumNDpy4WZ3Rtb9mPMBBNCwd0zlqtB+gU+bBJKeXaD5csdpwVkaN/eIiXB52sdE1ERtn/jOzpQSTR6pH56mlMAJmp5chWIZrhuvCuuw7OxyfpV2wSxX3oZB5Uy2FETA7veS/1HEm4+Xkp4XV3X12woV8nA9NsZKM3qrPCIrjm7eCxITW51sEpjg9eamLRFtO9coh8JsuQo59GwI5dG6w4ENV7sgA61oyETA5vQzmvjCyPl9LGIG3wim0ncAaLN617JOT37Aa9WkMAZPTi4kzQI30MUQ1WZXT9yasume1i5L9lVZ1SHE4HlLxsor0OLws2Tlyp+hndXoMQD6/yvqx5DCy9Iqe3ZrS0B6TBW9abFhjb5xfZ9IUIv+vusaJLeUfEQ8wsDCYQyt96ayXJj9Td1JAF1SyqtzMotoSka9iiLVK9m7FtRTKkSMjEBRKTELJ/93omiO6GI3a8EjM7ug8+sVd7+kUbE0v0d6e/uLJKMM3Dgq+90NH7/caovmEDjqy9yxlR3oxZXLGsUqUlxwKjSQ6Aiyy+Uwf3a2od/tW+K3DRtiv6vitrV/q7BxuDFdGsttjBFPeSTYOz2jZdCid2dRIEPii6x4euqdjQoePJIsl0qMEpTHE+eXif827FRfXaF0BblSl0685j+JrcG9eEmq2nMyUVoR4jDxDh97oIXlEeb5d0NF+8eqi76COR9PcZkd6UTunXzPSh3S+Ts8c1qHVoCVxeCznGt1Xpfd0REZhjiUK0Zyt1n70pvMP5nWtvujt0w90szi9G8KcIQbM6oOyCMiAzqEc9mUFBLhJJrl8jo2ftI4W4FLkEc7Nr/EcBj7TkHOj92Of0L+3pDdOSEwJn3eRezWa4QzZ1t2KReTvCBOH/RrqsEI3mkrJJqQ1yfRnjrKSEVveZYTZ6D35vV/QG/1HRylLpMfpKEwu1yoYe9JuxbUUXEEOQWa1SC9ZrJJE521yNtLc1ZatJeuiSDYOj/zYBBu9eT9aLE6/ZrRplFYj5+5WHG38FXlJkFkTe2SNBTZWSUYj6YK4MoFld29sGRW9k05vSW+aaCOBRJ5kER3gfdZhMRj//+Vrgjl525Ho+CiRk3sm2DVDZRFdpc9/SWEtgWgPcmlFt2akR1bQxtKZVHPRtCaGiWQjlxw/lvZQFX9I+pvF6WGaXaQcwLG2OqmoPvt7+DlijEeycbmC5+rq5IiPxZFdtStxVtuIoKyaFAk0lkZLJp0NX0htyYh22K48i9M3Q54E0/BgJOgyLnDtvsD5rNyIkxoR/HpX+b7BfqUGICsUUHD6B7I9KUa33NMNoy7PicSxfLtbJ/ta/FhNuqvzRYTXdZPC6Ox3QPROYmVxeufoWb0pu1VxQnYrhgA8UpXe6oLo1uqIcW2XewPevnGj0gpggazGqbpCxBfMe7jzXs9/8lznC6dsn8+2iXjSoehovh/QEULziiysbqBgZvm3pQtFgUwWJrKs3YoCBH1+6jz2ww2T4B0kznL6Qlgwgj4wYAD+PZ8ooYd40PEYPeg8bvOe8wkm02cBHUmF0X+0WJw+JGSU+ZyXahSGTu5WbBsPLeD0wJ8jf4fq0qbkkI2+hbD4W6yIsKNFsnBA97sKFhbe62i2JeiMr8y1CbxxfPJ5/4j5IlmA00mX4WJJbxzwDCGrFIEHcL39FyHBoiPnY8jbiSJdMJYbkt3cpLc7MEdR8XGdOGbkAZu5WKg++502Oti5GI2aejdqF72fZRPqDxo38kD0YL7ReZCC+dUskR4HpDTMzy936i+5X7ZzLKV5Cz5RMLk7MEMHh7VLLTOqNCRzscjauJWhzFw9MkRckBBk0DtZsp1eRnS9PraEmBS2FIF0IyJaHWMpdQYOI/0wCB3SX/lI6eCc3gWcFP2mBBbiwSRx1u9zS+iRbGxLIFnghL6hye2ATbbTtzUBn+F61wo/IIh9pVtbV/LK7N2BCeWc07uAkyBvXpfCgtWYOR+iSumR02b0HE23ukysTg/4zhBdPpd6vUMODw7klNbdgTH8mlUUtc8S6Yt0fslkZCJFm7Rnfi2btPqhZzZGlokshmFm3PznLQm59PwZhXatCJfeenVHv2ZFRO8tUZtRIzkx68goetk6UpFOVFSxDjPpaNbPGundL3kwRi2ZCMgwODwjpFvRSdZUHeiZqySX4JLMKkEPAZPTyxBEKSL8kx6XepNCgAns2O7AUH0XWPhdcC0FETA5fSv3nc7OIAVV2U9WChQ5uwNDjWdzFXOD7Fw1xKQ+GyIwxendp+Hq+EM8e0/ksLm7A3v12hs+ts1ZPQnxrM+EgNnpZchmdUJ1YxO2CqwQEE447WB3oJ4fpgASDZ2FkrXufE5a/+YiYHb6ljEfYFiHZaJVi4eAcGHiGtsd+MMjD90y4a3//z+ETIFnN6fwkFHZ9ovjv9SRvSIxRdZVqyNciOSkgGy39dMTNmaNLrHpPWkjPFI/mBBJLVMQmBrpkYWzP5eRarTvIz9ld2CXg4vydQLbRaXg9aRIj3w5O9GeT7812ncMIkxMuwM7VcGUKM/h//euLlm9nonAnEiP6D90kNszHNcyHwHSoWPhybeQWhZCYJbTyzh8rDrV4eeuC6m7v2yFJXMk9rLXj1ELm3mW06ObjMTE7Kw12sLq7id7YdekNTrXZeAVTHxwcnLCf4UiyvCl9aLIAKYhVvT8RvVU5xqpLlBMXwgvljf/0vG7FfM09+tNITz9RRa+fzx0kZ70hF+1u4NOYCoSwI9L2OLZ7CY0Vb6mxK1RqsMvZ3+Wjp1Pc278+n9hAErNt0BwhQAAAABJRU5ErkJggg==",
+      "text/latex": [
+       "$\\displaystyle \\left[ \\left( -5, \\  - \\frac{5}{2}\\right), \\  \\left( 4, \\  2\\right)\\right]$"
+      ],
+      "text/plain": [
+       "[(-5, -5/2), (4, 2)]"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "text/plain": [
+       "'Jacobian:'"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "data": {
+      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAIcAAAA2CAYAAADtcYraAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAGGklEQVR4Ae2d7VHcMBCGDyYFJKED0kGADi4dAKkgSQfJ8Av+MUkHJBUQ6ABSQQIdJKmAQAfkfYx1o7Nl38lfd5K1M0L2WpJ3V69317LPbBwfHz+fOOjk5OTBwU6sCC2guXZiYFO63qjcF8rnCG2QVKq2wEVh/sHDxYY8x29tfBB6rqv7hnFEOrxGKZWdsXg+6bktfT/lM7Sr+h/74t/mPO9Kfd+r08Ez755r1kGK4BK/qWAUjIOxRkE5MM5UvzEKaxuvfwNPpdUFT1gJmmSAB5UDlQ9S5DxoZfyFBwjoPSPZAS9CvogHbUXBg6OV9uF3nkqF3wJEMaHEYzwXv5UXTeAIGyCA4I9AgKdwURE0rjaVvOBzjkrNRnBAoDioUJPEfKLjjZNS+ifPgRUiIgECYBBOzB1MY+0SOBqbbm07koheCiRf2kq48rAiJYiLP1R84iN3J61cZlvDrWN/2eRMcpGDVIUbL7HXARwkUzteUqfGJQsIECxcvVQ9W/MoNfJkpLDiabB1bC5A7EuuV7bH0PY2pY28CRxtrLcGfQUAEtA91cUEFMCwatyYVh5WGkvu7riVs1+qJlxFTblnIAG91jb5hk1T8VolpVGAQ0YwS8XT3DoX4v3R9pXqr7bFItu+kj6EDvKNIrVO2GMBRyfZedG6674v4L/qU8aUc/Rp3cDHTuAIfAL7FD+Bo0/rBj52AkfgE9in+AkcfVo38LETOAKfwD7FD/ZWVrdxPKg7VOHl6NKzGfHMG/R3arOl/eIKYp92jWLsID2HJpolYxZ+AAhljnScY3cAQoVVQhbDWDBK5GGBUMFxm086q6AuwmvM3rxWW7ZZTt52NU48twWCBIdblSduDgC8SRE4D+LhcRItaYFgc44a/UphJm/LE0oeyI2SdNFwYfAMaukffMUIjqrJBxhVwKnqEzRfgEDfxj/4ihEchA8XYahiqHG1i4YncGCL7KGktj9q2yusxphzAACM4ko+Wz/G1rijoejAkc/cqeqpmUVdNWzzRvaoPIfRv2kdZFjRJOMV+I0ok867kty6sq6RvflEDc/wdYz3Ht6pJPKwQKjgwAPUrngKGLXHPWw02qaxhpXRTmiXitd6Dl19ZPjpB0ddWtxjrFXbfxE4yPpLD7U89MuaSslH3z6htJduG33JqrE7sX9T+WrB0XTQYr8+DVg8V9f7kn2qMfGgJMEktrzZPntuo/1oaRBwBG49wPACHVQDEj6m1pu34DzrQitJSGXkfRWuxFrKJ6O2TRcHdR5WD6vIDqvIjKsPkewffC0l/+Dg0ETwrgXf8apdkMqB8Ve115JvldaMp/JehU9rzpF4rIsUfzGWtRHflpP1lKB+IyP58Xw8cDM/fMr2xTP7c7awdwYNKxIIt8xX7pYx8JHamlhfu+yt8fii3tyH02wldQyAkTtAjOkijMbCmXN9RHzA0+oTjq6T9s2T3MvY2inG0J6DK895hdrSSSEm8lwFF74w/Cxqo/EWvRw0URvzQlAJPDqG3BkwtF0XgtQsHhoaHLyNtUym/1rt8Ba49L0BzZ0BxD6f5CAMAYh7bXNLXumh7H4xbA8WVmRYPACeoJbU7qOK+XX4LzXmw7ND0U+diI+fXJoTShY7ITXsUdRDeg7AwdtYlaSJIJzYnoVPb9NvKFo2jA0lz0rPM5jnkJbE8kWeg+SJO5m3uVUABncZlKyvanKWImB2xS++XU6C6fv5hVrw5jKNphoSHOQPAMRJmkjCyVw81z4g4As11NkdS7ENg4nHTw+6+BYW57FvXRl+tLQ5lOaaPCYX45dIxwgnJa8ivpkoZ7/SQO0ZnIdQlkgWMJ6DF2bmrmrtlyarA4vxGUTOZSZ9om3CxKEKucYsDIiPPNxCQkfaJ9zY+cjTkSX+qh+TXvlykDUE3mfOe1nHot3MbW3rl2HBgIMJstcfmKQ+jMRkM+5soUmCsV86l/iA03lMfC/SWIBxdk5XZ2OgvK2rScw8VlDx3jZdTx4fHwct+uc/ZyrbXZ6XMduO14dcbWVadf8NBBiadHXiQU5V9xG6vNWRHDxnIOQ1ClveJwykw+Yq5NQk1Lr4Fcj0PQGjbPX/J2hbH4mZLlwAAAAASUVORK5CYII=",
+      "text/latex": [
+       "$\\displaystyle \\left[\\begin{matrix}-1 & 2\\\\- \\frac{10}{\\left(A + 1\\right)^{2}} & -1\\end{matrix}\\right]$"
+      ],
+      "text/plain": [
+       "⎡   -1     2 ⎤\n",
+       "⎢            ⎥\n",
+       "⎢  -10       ⎥\n",
+       "⎢────────  -1⎥\n",
+       "⎢       2    ⎥\n",
+       "⎣(A + 1)     ⎦"
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    },
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Equilibrium point (-5, -5/2) is classified as: Stable Spiral\n",
+      "Equilibrium point (4, 2) is classified as: Stable Spiral\n"
+     ]
+    }
+   ],
+   "source": [
+    "A, B = sp.symbols('A B')\n",
+    "\n",
+    "L = 10\n",
+    "f = -A + 2 * B\n",
+    "g = -B + L / (A+1)\n",
+    "F = sp.Matrix([f, g])\n",
+    "X = sp.Matrix([A, B])\n",
+    "\n",
+    "eq_points = equilibrium_points(F, X)\n",
+    "display(\"Equilibrium point:\", eq_points)\n",
+    "Jacobian = find_jacobian(F, X)\n",
+    "\n",
+    "display(\"Jacobian:\", Jacobian)\n",
+    "for X_eq in eq_points:\n",
+    "    classification = classify_equilibrium(Jacobian, X_eq, X)\n",
+    "    print(f\"Equilibrium point {X_eq} is classified as: {classification}\")\n"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,