Laplace_equation en.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "22261ebd-fe17-445e-ae5f-33fff77d2b28",
   "metadata": {},
   "source": [
    "# Laplace Equation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cf0611db-f2b5-4cd1-94a9-61f15111e170",
   "metadata": {},
   "source": [
    "The Laplace's equation is a second-order partial differential equation named after [Pierre-Simon Laplace](https://en.wikipedia.org/wiki/Pierre-Simon_Laplace), who first studied its properties. The homogenous and its inhomogeneous counterpart, the [Poisson equation](https://en.wikipedia.org/wiki/Poisson%27s_equation), are governing the mechanical and thermal equilibria of bodies, as well as ﬂuid-mechanical and electromagnetic potentials. The Laplace equation is a homogeneous linear diﬀerential equation.\n",
    "\n",
    "As we did with the diffusion equation, we will use the method of separation of variables to solve the homogeneous Laplace equation in a particular boundary-value problem. In the cartesian coordinates, the Laplace equation in 3-D is given by:\n",
    "\n",
    "$$ \\nabla^2 V \\space = \\space \\frac{\\partial^2 \\space V}{\\partial x^2} + \\frac{\\partial^2 \\space V}{\\partial y^2} + \\frac{\\partial^2 \\space V}{\\partial z^2}$$\n",
    "\n",
    "The theory of potential allows us to solve for V from given boundary conditions. After we obtain the potential V, we can ﬁnd the electric ﬁeld $\\vec{E}$ by calculating the gradient of the potential,\n",
    "\n",
    "$$\\vec{E} \\, = \\, - \\nabla \\, V $$\n",
    "\n",
    "Here is an illustrative example with a two-dimensional situation taken in Griffiths, D. J.,*Introduction to Electrodynamics*, Fourth Edition, Cambridge University Press, 2017, example 3.4, page 136\n",
    "\n",
    "<hr style=\"border:1px solid gray\">\n",
    "\n",
    "Four long metallic plates form a shaft along the *z* direction; two of them are grounded at $y=0$ and $y=a$, and the other two are maintained at a constant potential $V_0$ at $x \\pm b$. See the following figure.\n",
    "\n",
    "![metallic plate](Example-3.4.png)\n",
    "\n",
    "So the boudary conditions are:\n",
    "\n",
    "![Boundary conditions](Boundary_conditions.png)\n",
    "\n",
    "<hr style=\"border:1px solid gray\">\n",
    "\n",
    "Our task is to find the potential inside the volume containt by the four plates. To solve this problem, we assume that *V* is independant of *z*. The method we will use is the separation of variables, i.e. we use the ansatz\n",
    "\n",
    "$$  V = X(x) \\space Y(y)$$\n",
    "\n",
    "where $X$ and $Y$ are functions of a single variable *x* and *y*, respectively. Start a new Jupyter notebook and define the symbols for our problem:\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "d0c06d2b-28e8-484c-b232-def50774b24b",
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import *\n",
    "\n",
    "x, y = symbols(\"x, y\", real=True)\n",
    "\n",
    "V, X, Y = symbols(\"V, X, Y\", cls=Function)\n",
    "\n",
    "V = V(x, y)\n",
    "X = X(x)\n",
    "Y = Y(y)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5d63ba7c-4b39-425a-a11c-fba39c0d1cc8",
   "metadata": {},
   "source": [
    " We also define a convenience function `apply` which allows to apply a given operation to both sides of an equation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "a4f418e8-ee8b-4a20-a0eb-326dc2a3b63b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def apply(eq, func, *args):\n",
    "    return Eq(\n",
    "        func(eq.lhs, *args), \n",
    "        func(eq.rhs, *args)\n",
    "    )"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5998242e-d57b-42a8-bad3-9f9dad2c954c",
   "metadata": {},
   "source": [
    "Now we can formulate the partial equation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "2bba8c2b-3719-479b-9e0d-9494bfdba8e1",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ1 = Eq(diff(V,(x,2)) + diff(V, (y,2)), 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "15963769-608f-4640-a822-f54e52899665",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\partial^{2}}{\\partial x^{2}} V{\\left(x,y \\right)} + \\frac{\\partial^{2}}{\\partial y^{2}} V{\\left(x,y \\right)} = 0$"
      ],
      "text/plain": [
       "Eq(Derivative(V(x, y), (x, 2)) + Derivative(V(x, y), (y, 2)), 0)"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ1"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac91bb0d-85b8-4b62-8dc5-c37a79f0a511",
   "metadata": {},
   "source": [
    "Now, we substitute our separation ansatz:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "d0bdc42a-f619-4bfd-9aa8-339ab1d8a61e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle X{\\left(x \\right)} \\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)} + Y{\\left(y \\right)} \\frac{d^{2}}{d x^{2}} X{\\left(x \\right)} = 0$"
      ],
      "text/plain": [
       "Eq(X(x)*Derivative(Y(y), (y, 2)) + Y(y)*Derivative(X(x), (x, 2)), 0)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ1.subs(V, X * Y).doit()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a073359d-89cb-4601-aaa4-303ca5e8b6e9",
   "metadata": {},
   "source": [
    "The derivative operator has been pulled through the factors that were not affected by the derivative. And since the remaining derivative only acts on a function of a single variable, *sympy* has replaced the partial derivative with a total derivative.\n",
    "\n",
    "The next step in the separation method is to divide by the product *X(x)Y(y)*:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "878bf5ef-a5ae-4d4e-ba7d-eef16c130056",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{X{\\left(x \\right)} \\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)} + Y{\\left(y \\right)} \\frac{d^{2}}{d x^{2}} X{\\left(x \\right)}}{X{\\left(x \\right)} Y{\\left(y \\right)}} = 0$"
      ],
      "text/plain": [
       "Eq((X(x)*Derivative(Y(y), (y, 2)) + Y(y)*Derivative(X(x), (x, 2)))/(X(x)*Y(y)), 0)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "apply(_, Mul, 1 / (X * Y))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "f03deef7-2ecf-4546-958f-fd44272db304",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)}}{Y{\\left(y \\right)}} + \\frac{\\frac{d^{2}}{d x^{2}} X{\\left(x \\right)}}{X{\\left(x \\right)}} = 0$"
      ],
      "text/plain": [
       "Eq(Derivative(Y(y), (y, 2))/Y(y) + Derivative(X(x), (x, 2))/X(x), 0)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "expand(_)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4dbfa441-20ee-459a-9f64-8e753a2e0cc9",
   "metadata": {
    "tags": []
   },
   "source": [
    "We have now to put the two terms on different sides of the equation, adding the negative of the first (`args[1]`) term of the left-hand side:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "0ff0354b-7b8d-4047-b780-fb70fcf7819d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\frac{d^{2}}{d x^{2}} X{\\left(x \\right)}}{X{\\left(x \\right)}} = - \\frac{\\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)}}{Y{\\left(y \\right)}}$"
      ],
      "text/plain": [
       "Eq(Derivative(X(x), (x, 2))/X(x), -Derivative(Y(y), (y, 2))/Y(y))"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "apply(_, Add, - _.lhs.args[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "0ff196db-39d9-4f56-9dbb-2a3aca1e2b08",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ2 = _"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ca7a41f2-4001-4719-930e-04aaa0ac6d5e",
   "metadata": {},
   "source": [
    "The key observation here is that the left side depends only on *x*, whereas the right side depends only on *y*. So, when you change *x*, only the left side changes but not the right side. Conversely, if you change *y*, the right side changes, but not the left side. The only way how this can be true is that both sides are equal to a constant. Let’s call that constant −𝑘²:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "c2c9a870-72e4-44eb-a66d-daf28523dd99",
   "metadata": {},
   "outputs": [],
   "source": [
    "k = symbols(\"k\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "6ee31bdb-5dc9-4b83-83e9-adfd51b46305",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle - \\frac{\\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)}}{Y{\\left(y \\right)}} = k^{2}$"
      ],
      "text/plain": [
       "Eq(-Derivative(Y(y), (y, 2))/Y(y), k**2)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(EQ2.rhs, k**2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e9537238-48f6-4c47-9f09-4d61a2bdce95",
   "metadata": {},
   "source": [
    "Lets multiply the result by *X(x)*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "5cb9fb10-59b3-4f1d-ada8-d643f05f0909",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle - \\frac{d^{2}}{d y^{2}} Y{\\left(y \\right)} = k^{2} Y{\\left(y \\right)}$"
      ],
      "text/plain": [
       "Eq(-Derivative(Y(y), (y, 2)), k**2*Y(y))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "apply(_, Mul, Y)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2b2aa3d1-e44e-4a75-9a61-92bb4f72f87f",
   "metadata": {},
   "source": [
    "This is the ordinary differential equation of a harmonic oscillator. It can be easily solve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "ad52e700-f293-43ed-afcc-6202160df8da",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ3 = _"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "a99a5e43-8ae3-48c1-aa66-1efc6ab9c08b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle Y{\\left(y \\right)} = C_{1} e^{- i k y} + C_{2} e^{i k y}$"
      ],
      "text/plain": [
       "Eq(Y(y), C1*exp(-I*k*y) + C2*exp(I*k*y))"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol_Y = dsolve(EQ3, Y)\n",
    "sol_Y"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "478e5faf-327d-48b2-8ee0-5925dc75de81",
   "metadata": {},
   "source": [
    "Of course, we could have written directly the solution for $Y(y)$ (See: diffusion_equation en.ipynb)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "2a497fbf-2c56-4fce-b33f-e7bbb4aa736d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(A, B)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var(\"A, B\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "f1bd8627-b512-49e9-9828-5dc710810c45",
   "metadata": {},
   "outputs": [],
   "source": [
    "sol_Y = A*sin(k*y) + B * cos(k*y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "bf6f410e-c68e-4a1b-985e-4238947a9b42",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle A \\sin{\\left(k y \\right)} + B \\cos{\\left(k y \\right)}$"
      ],
      "text/plain": [
       "A*sin(k*y) + B*cos(k*y)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol_Y"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6f3cedbf-924c-4ac9-93ea-c419b6c786f4",
   "metadata": {},
   "source": [
    "The boundary condition that $V = 0$ when $ y = 0$ give us"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "1b03c396-4b26-459c-9b71-bf7910b0e019",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle A \\sin{\\left(k y \\right)}$"
      ],
      "text/plain": [
       "A*sin(k*y)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "_.subs(B,0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "ca5bffe7-a1bc-4d37-989a-c2bd65fd1002",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ4 = _"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0ae193c2-f3e7-46d4-88f0-547927dfeeac",
   "metadata": {},
   "source": [
    "The other boundary condition gives a constraint on *k*:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "4adc4e5b-4d1f-4ffa-a82d-ba984d6e6d83",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle a$"
      ],
      "text/plain": [
       "a"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var(\"a\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "eee25790-b771-4778-8985-5ee9051cb603",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle A \\sin{\\left(a k \\right)} = 0$"
      ],
      "text/plain": [
       "Eq(A*sin(a*k), 0)"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(EQ4.subs(y,a),0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f39189f1-6b95-4542-8de5-e0f277b2f365",
   "metadata": {},
   "source": [
    "Since the constant *A* cannot be 0, we must instead ask that $a k = n \\pi$. We know that this is 0 for all integer mulitples of $2 \\pi$. So lets define that variable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "bb83f15c-3f1e-4c2d-adb7-c82c58653fcc",
   "metadata": {},
   "outputs": [],
   "source": [
    "n = Symbol(\"n\", integer=True, positive=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "6e70f273-bd86-4bfb-a782-048979773cbf",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle k = \\frac{\\pi n}{a}$"
      ],
      "text/plain": [
       "Eq(k, pi*n/a)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(k,n*pi/a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "966989ad-d841-4e40-9a9c-71bc648a6743",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ_k = _"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "efee77dc-5ef8-4ea6-a97a-627e82ef78ac",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle A \\sin{\\left(\\frac{\\pi n y}{a} \\right)}$"
      ],
      "text/plain": [
       "A*sin(pi*n*y/a)"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ4.subs(k,EQ_k.rhs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "45a783a4-45bd-419b-b799-ac94d110ee4c",
   "metadata": {},
   "outputs": [],
   "source": [
    "sol_Y = _"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "436f84b1-9219-43d6-bf00-6e88f28f57f0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle A \\sin{\\left(\\frac{\\pi n y}{a} \\right)}$"
      ],
      "text/plain": [
       "A*sin(pi*n*y/a)"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol_Y"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4585e685-e334-4f3f-b513-163eaa130369",
   "metadata": {},
   "source": [
    "Now, lets solve for $X(x)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "030be343-3e9f-4ff3-9781-b798e45df9fa",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\frac{d^{2}}{d x^{2}} X{\\left(x \\right)}}{X{\\left(x \\right)}} = k^{2}$"
      ],
      "text/plain": [
       "Eq(Derivative(X(x), (x, 2))/X(x), k**2)"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(EQ2.lhs, k**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "6689aa16-4f46-4373-9b71-59337bc26c26",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ5 = _"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "id": "de7b5dc8-1290-450b-8f73-12b776b0da33",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle X{\\left(x \\right)} = C_{1} e^{- k x} + C_{2} e^{k x}$"
      ],
      "text/plain": [
       "Eq(X(x), C1*exp(-k*x) + C2*exp(k*x))"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol_X = dsolve(EQ5, X)\n",
    "sol_X"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2fd34657-4ad6-4937-9894-b85d281fdd7c",
   "metadata": {},
   "source": [
    "Because the boundary conditions on x is restricted on the range $x  \\in [-b,b]$ , the exponentials are perfectly acceptable.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "210038c8-acfe-4d87-a8b4-ec3536a80881",
   "metadata": {},
   "source": [
    "Now, since the situation is symmetric with respect to x, we get that $X(-b)=X(b)$, it follow that $C_1 = C_2$. To get around this, we are going to be using the fact that\n",
    "\n",
    "$$\n",
    "e^{k x}+e^{-k x}=2 \\cosh(k x)\n",
    "$$\n",
    "\n",
    "and will be absorbing the 2 and $C_1$ with the new variable $C$. Now we can write a solution for $X(x)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "id": "b0462755-73e1-4fb8-8ba1-43e8fc1c02c3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle C$"
      ],
      "text/plain": [
       "C"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var(\"C\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "id": "415cdd91-aab7-4828-a547-5d750e6e0864",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle C \\cosh{\\left(k x \\right)}$"
      ],
      "text/plain": [
       "C*cosh(k*x)"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sol_X = C*cosh(k*x)\n",
    "sol_X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "fcad1c90-a0dd-4d19-83a5-ae0a9687705a",
   "metadata": {},
   "outputs": [],
   "source": [
    "EQ6 = _"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "id": "46a1e9d6-6eeb-48b7-80f1-6b7bba63573e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle C \\cosh{\\left(\\frac{\\pi n x}{a} \\right)}$"
      ],
      "text/plain": [
       "C*cosh(pi*n*x/a)"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ6.subs(k,EQ_k.rhs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "id": "b6f9d897-1614-4cd9-ad93-aaf912747407",
   "metadata": {},
   "outputs": [],
   "source": [
    "sol_X = _"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa871a6a-fb2e-4e27-9090-09494c3d266d",
   "metadata": {},
   "source": [
    "Now we can assemble a solution $V$ from our base solutions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "id": "2d39c520-3ccb-47d7-92ee-89a1f7a1bee5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle V{\\left(x,y \\right)} = \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi n x}{a} \\right)} {a}_{n}$"
      ],
      "text/plain": [
       "Eq(V(x, y), Sum(sin(pi*n*y/a)*cosh(pi*n*x/a)*a[n], (n, 1, oo)))"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a = IndexedBase(\"a\")\n",
    "\n",
    "sol_V = Eq( V, Sum( sol_X * sol_Y * a[n], (n, 1, oo)) ).subs({A:1,C:1})\n",
    "sol_V"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "07592c2e-128b-41a7-b011-2246c353831c",
   "metadata": {},
   "source": [
    "Now, we have to take into account to the initial condition on$x$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "id": "89b977fa-121a-4229-97f5-6cc854933680",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle V{\\left(b,y \\right)} = \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}$"
      ],
      "text/plain": [
       "Eq(V(b, y), Sum(sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo)))"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var(\"b\")\n",
    "_.subs(x,b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "id": "6144d7af-ef74-4ba2-9a40-e03019380993",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}$"
      ],
      "text/plain": [
       "Sum(sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo))"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ7 = _.rhs\n",
    "EQ7"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "id": "26609068-b63c-4ebe-a96b-04f9cbca1eb5",
   "metadata": {},
   "outputs": [],
   "source": [
    "Vo = symbols(\"Vo\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "id": "0dfac985-0918-4cf7-aaff-e2a27b6d261b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n} = Vo$"
      ],
      "text/plain": [
       "Eq(Sum(sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo)), Vo)"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(EQ7,Vo)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93f65d6b-38e6-42f3-8301-e77d3098dc65",
   "metadata": {},
   "source": [
    "This expression is simply a Fourier series. We can solve for $a_n$ by using the fact that the sines are *orthogonal* functions. That mean that their scalar product satisfies\n",
    "\n",
    "$$\n",
    "\\int_0^{2 \\pi} \\sin(  n x) \\sin(m x)  dx =\\begin{cases}0 & n  \\neq  m\\\\\\pi & otherwise\\end{cases} \n",
    "$$\n",
    "\n",
    "So we multiply both sides by a suitable sine and integrate:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "id": "ecbe43fb-9706-4c18-a97d-fc28145a6c4c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n} = Vo$"
      ],
      "text/plain": [
       "Eq(Sum(sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo)), Vo)"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ8 = _\n",
    "EQ8"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "id": "e6d72e3d-67ac-468c-a2ba-3b27e97bd218",
   "metadata": {},
   "outputs": [],
   "source": [
    "m = Symbol(\"m\", integer=True, positive=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "751fc78c-dc17-4bba-bb15-785ce6207066",
   "metadata": {},
   "source": [
    "and now we multiply both sides of EQ8 by $\\sin  \\big( \\frac{ \\pi m y}{a}\\big) $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "id": "3cc0e37c-9024-42d6-b6a8-cc99489c1425",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi m y}{a} \\right)} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}$"
      ],
      "text/plain": [
       "Sum(sin(pi*m*y/a)*sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo))"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Sum(a[n]*sin(pi*n*y/a)*sin(pi*m*y/a)*cosh(pi*b*n/a),(n,1,oo))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "id": "a781ada9-1527-44de-9eb8-7d0b6b173b1f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi m y}{a} \\right)} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}$"
      ],
      "text/plain": [
       "Sum(sin(pi*m*y/a)*sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo))"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ9 = _\n",
    "EQ9"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "id": "56d4da03-f274-4452-84dd-80be1eb0cb2e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle Vo \\sin{\\left(\\frac{\\pi m y}{a} \\right)}$"
      ],
      "text/plain": [
       "Vo*sin(pi*m*y/a)"
      ]
     },
     "execution_count": 45,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Vo*sin(pi*m*y/a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "6f3257eb-5c55-4f13-8e88-9cab49fa174c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle Vo \\sin{\\left(\\frac{\\pi m y}{a} \\right)}$"
      ],
      "text/plain": [
       "Vo*sin(pi*m*y/a)"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "EQ10 = _\n",
    "EQ10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "id": "9b446013-65c1-4734-b009-315b6adeb4f7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\int\\limits_{0}^{a} \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi m y}{a} \\right)} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}\\, dy$"
      ],
      "text/plain": [
       "Integral(Sum(sin(pi*m*y/a)*sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo)), (y, 0, a))"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Integral(EQ9,(y,0,a))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "id": "282614b9-f709-447c-ab0a-bdc47d7f6465",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\int\\limits_{0}^{a} \\sum_{n=1}^{\\infty} \\sin{\\left(\\frac{\\pi m y}{a} \\right)} \\sin{\\left(\\frac{\\pi n y}{a} \\right)} \\cosh{\\left(\\frac{\\pi b n}{a} \\right)} {a}_{n}\\, dy$"
      ],
      "text/plain": [
       "Integral(Sum(sin(pi*m*y/a)*sin(pi*n*y/a)*cosh(pi*b*n/a)*a[n], (n, 1, oo)), (y, 0, a))"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "_.doit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1b5f70d2-3382-4678-8cb3-8dc9e5d30df7",
   "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.9.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}