Chapter-15.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib widget\n",
    "from sympy import *\n",
    "init_printing(use_latex=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 15 - The Equation class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from algebra_with_sympy import *\n",
    "# by default, when an equation is shown on the screen, it also shows\n",
    "# a label with a unique number. Usually, I hide them:\n",
    "algwsym_config.output.label = False\n",
    "# by default, the algebra_with_sympy's solve returns objects of type\n",
    "# FiniteSet. I don't like that behavior: I prefer lists of solutions.\n",
    "algwsym_config.output.solve_to_list = True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a, b, c, x = symbols(\"a, b, c, x\")\n",
    "eq1 = Equation(a + b, c * x)\n",
    "eq2 = Equation(x**2 - 4 * x, a + b)\n",
    "display(eq1, eq2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "display(eq1.lhs, eq1.rhs, eq1.reversed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "res = solve(eq1, x)\n",
    "res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1 + 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1 / eq2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1.applyrhs(lambda t: t**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1.subs(c, 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eq1.subs(eq2.reversed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 15.1 - Electric Circuit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "L, R, C1, C2, t = symbols(\"L, R, C1, C2, t\")\n",
    "iC1, iC2, VL, Vi = symbols(\"i_{C_{1}}, i_{C_{2}}, V_L, V_i\")\n",
    "VC1, VC2, iL = [f(t) for f in symbols(\"V_{C_{1}}, V_{C_{2}}, i_L\", cls=Function)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "iC2_eq = Equation(iC2, (VL - VC2) / R)\n",
    "iC1_eq = Equation(iC1, iL + iC2)\n",
    "VL_eq = Equation(VL, Vi - VC1)\n",
    "display(iC1_eq, iC2_eq, VL_eq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "e1 = Equation(C1 * VC1.diff(t), iC1)\n",
    "e2 = Equation(L * iL.diff(t), VL)\n",
    "e3 = Equation(C2 * VC2.diff(t), iC2)\n",
    "display(e1, e2, e3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "e1 = (e1.subs(iC1_eq).subs(iC2_eq).subs(VL_eq) / C1).expand()\n",
    "e2 = (e2.subs(VL_eq) / L).expand()\n",
    "e3 = (e3.subs(iC2_eq).subs(VL_eq) / C2).expand()\n",
    "display(e1, e2, e3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A, b = linear_eq_to_matrix([e.rhs for e in [e1, e2, e3]], [VC1, VC2, iL])\n",
    "A, b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "final = Equation(\n",
    "    Matrix([e.lhs for e in [e1, e2, e3]]),\n",
    "    MatMul(A, b, evaluate=False))\n",
    "final"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 15.2 - Temperature Distribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "t, k, x, m, Cp, rho, alpha = symbols(r\"t, k, x, m, C_p, rho, alpha\", positive=True)\n",
    "Dt, Dx, Dy, Dz, DV = [symbols(r\"{\\Delta}\" + s, positive=True) for s in \"t,x,y,z,V\".split(\",\")]\n",
    "Q, q, P, T = [e(t, x) for e in symbols(\"Q, q, P, T\", cls=Function)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Mass = Equation(rho, m / DV)\n",
    "Veq = Equation(DV, Dx * Dy * Dz)\n",
    "Qeq = Equation(Q, m * Cp * T)\n",
    "fourier = Equation(q, -k * T.diff(x))\n",
    "td = Equation(alpha, k / (Cp * rho))\n",
    "display(Mass, Veq, Qeq, fourier, td)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# change in heat content\n",
    "hc_change = Q - Q.subs(t, t + Dt)\n",
    "# heat in: Dy * Dz is the surface area through which q flows\n",
    "hi = q * Dt * Dy * Dz\n",
    "# heat out: Dy * Dz is the surface area through which q flows\n",
    "ho = q.subs(x, x + Dx) * Dt * Dy * Dz\n",
    "# heat generated\n",
    "hg = P * Dt * Dx * Dy * Dz\n",
    "# heat balance equation\n",
    "hb = Equation(hc_change, hi - ho + hg)\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Veq = Equation(DV, Dx * Dy * Dz)\n",
    "hb = hb.applyrhs(lambda e: ((e * Dx).expand().subs(*Veq.reversed.args) / Dx)).expand()\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = (hb / Dt).expand()\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = hb.applylhs(lambda e: Limit(e.collect(1 / Dt), Dt, 0))\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = hb.applylhs(lambda e: e.subs(e, Q.diff(t)))\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "qxdx = Equation(q.subs(x, x + Dx), q.subs(x, x + Dx).series(Dx, 0, 2).removeO().doit())\n",
    "qxdx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = hb.subs(qxdx).expand()\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = hb.subs(fourier).subs(Qeq).doit()\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = (hb / DV / k).expand().subs(Mass.reversed)\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hb = hb.subs(1 / td.reversed)\n",
    "hb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.12"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}