From 31369b24398d2f6798c45d50398364c1d4df12b2 Mon Sep 17 00:00:00 2001 From: Ziaeemehr Date: Sat, 10 Aug 2024 19:33:49 +0200 Subject: [PATCH] chapter 4 updated --- examples/chap_04.ipynb | 124 +++++++++++++---------------------------- 1 file changed, 38 insertions(+), 86 deletions(-) diff --git a/examples/chap_04.ipynb b/examples/chap_04.ipynb index 5e05594..d29cc48 100644 --- a/examples/chap_04.ipynb +++ b/examples/chap_04.ipynb @@ -16,13 +16,14 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import sympy as sp\n", "import numpy as np\n", - "from scipy.linalg import eig" + "from scipy.linalg import eig\n", + "from IPython.display import display, Math" ] }, { @@ -33,87 +34,47 @@ ] }, { - "cell_type": "code", - "execution_count": 2, + "cell_type": "markdown", "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "array([-3.5+14.72243186j, -3.5-14.72243186j, -8. +0.j ])" - ] - }, - "execution_count": 2, - "metadata": {}, - "output_type": "execute_result" - } - ], "source": [ - "A = np.array([[-5, -10, 7], [7, -5, -10], [-10, 7, -5]])\n", - "B = np.array([1,1,1])\n", - "X0 = np.array([0, 0, 0])\n", - "\n", - "eigenvalues, eigenvectors = eig(A)\n", - "eigenvalues" + "$$\n", + "\\frac{d}{dt}\n", + "\\begin{pmatrix}\n", + "E_1 \\\\\n", + "E_2 \\\\\n", + "E_3\n", + "\\end{pmatrix} =\n", + "\\begin{pmatrix}\n", + "-5 & -10 & 7 \\\\\n", + " 7& -5& -10 \\\\\n", + " -10& 7 & -5 \n", + "\\end{pmatrix} \n", + "\\begin{pmatrix}\n", + "E_1 \\\\\n", + "E_2 \\\\\n", + "E_3\n", + "\\end{pmatrix}\n", + "$$" ] }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 2, "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "x1(t) = exp(-8*t) + 4*sqrt(3)*exp(-7*t/2)*sin(17*sqrt(3)*t/2)\n", - "x2(t) = exp(-8*t) - 2*sqrt(3)*exp(-7*t/2)*sin(17*sqrt(3)*t/2) - 6*exp(-7*t/2)*cos(17*sqrt(3)*t/2)\n", - "x3(t) = exp(-8*t) - 2*sqrt(3)*exp(-7*t/2)*sin(17*sqrt(3)*t/2) + 6*exp(-7*t/2)*cos(17*sqrt(3)*t/2)\n" - ] - } - ], + "outputs": [], "source": [ - "# Define the symbolic variables\n", - "t = sp.Symbol('t')\n", - "x1, x2, x3 = sp.Function('x1')(t), sp.Function('x2')(t), sp.Function('x3')(t)\n", - "\n", - "# Define the matrix A\n", - "A = sp.Matrix([[-5, -10, 7],\n", - " [7, -5, -10],\n", - " [-10, 7, -5]])\n", - "\n", - "# Define the vector X\n", - "X = sp.Matrix([x1, x2, x3])\n", - "\n", - "# Define the system of differential equations\n", - "dX_dt = A * X\n", - "\n", - "# Create the system of differential equations\n", - "system = [sp.Eq(sp.diff(x1, t), dX_dt[0]),\n", - " sp.Eq(sp.diff(x2, t), dX_dt[1]),\n", - " sp.Eq(sp.diff(x3, t), dX_dt[2])]\n", - "\n", - "# Solve the system of differential equations\n", - "solution = sp.dsolve(system)\n", - "\n", - "# Define the initial conditions\n", - "X0 = sp.Matrix([1, -5, 7])\n", - "\n", - "# Apply the initial conditions to find the constants\n", - "t0 = 0\n", - "constants = sp.solve([sol.rhs.subs(t, t0) - x0 for sol, x0 in zip(solution, X0)])\n", - "\n", - "# Substitute the constants back into the solution\n", - "final_solution = [sol.rhs.subs(constants) for sol in solution]\n", + "A = np.array([[-5, -10, 7], [7, -5, -10], [-10, 7, -5]])\n", + "B = np.array([0,0,0])\n", + "X0 = np.array([1, -5, 7])\n", + "t_range = np.linspace(0, 2, 1000)\n", "\n", - "# Print the final solution\n", - "for i, sol in enumerate(final_solution, 1):\n", - " print(f\"x{i}(t) = {sol}\")" + "# eigenvalues, eigenvectors = eig(A)\n", + "# eigenvalues" ] }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 3, "metadata": {}, "outputs": [ { @@ -154,22 +115,24 @@ } ], "source": [ - "from IPython.display import display, Math\n", + "from spikes.solver import solve_system_of_equations\n", + "\n", + "final_solution, x_values = solve_system_of_equations(A, B, X0, t_range)\n", "\n", "for i, sol in enumerate(final_solution, 1):\n", - " display(Math(f\"x_{i}(t) = {sp.latex(sol)}\"))" + " display(Math(f'x_{i}(t) = {sp.latex(sol)}'))\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "- Evaluating over a given time interval:" + "- plot the evaluation of analytical solution over a given time interval:" ] }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 4, "metadata": {}, "outputs": [ { @@ -188,16 +151,6 @@ "from sympy import lambdify\n", "import matplotlib.pyplot as plt\n", "\n", - "# Create a time range\n", - "t_range = np.linspace(0, 2, 1000) # From 0 to 10 seconds, with 1000 points\n", - "\n", - "# Convert symbolic solutions to numerical functions\n", - "x_functions = []\n", - "x_values = []\n", - "for sol in final_solution:\n", - " x_functions.append(lambdify(t, sol, 'numpy'))\n", - " x_values.append(x_functions[-1](t_range))\n", - "\n", "plt.figure(figsize=(6, 4))\n", "for i in range(3):\n", " plt.plot(t_range, x_values[i], label=f\"x{i+1}(t)\")\n", @@ -218,7 +171,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 5, "metadata": {}, "outputs": [ { @@ -249,7 +202,6 @@ "t = np.linspace(0, 2, 1000)\n", "solution_numerical = odeint(system, X0, t, args=(A,))\n", "\n", - "\n", "plt.figure(figsize=(6, 4))\n", "for i in range(3):\n", " plt.plot(t_range, x_values[i], label=f\"x{i+1}(t)\", alpha=0.5)\n",