Skip to content

Commit

Permalink
chapter 4 updated
Browse files Browse the repository at this point in the history
  • Loading branch information
Ziaeemehr committed Aug 10, 2024
1 parent 5f3bc26 commit 31369b2
Showing 1 changed file with 38 additions and 86 deletions.
124 changes: 38 additions & 86 deletions examples/chap_04.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -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"
]
},
{
Expand All @@ -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": [
{
Expand Down Expand Up @@ -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": [
{
Expand All @@ -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",
Expand All @@ -218,7 +171,7 @@
},
{
"cell_type": "code",
"execution_count": 6,
"execution_count": 5,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -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",
Expand Down

0 comments on commit 31369b2

Please sign in to comment.