{ "cells": [ { "cell_type": "markdown", "id": "4bc9f590-888f-42b9-81f7-1096612adde8", "metadata": {}, "source": [ "## The heat equation in a 1D Cartesian system (Dirichlet BC)\n", "\n", "We compare the FVM solution by PyFVTool to the analytic solutions, for two different types of boundary conditions." ] }, { "cell_type": "code", "execution_count": 1, "id": "106afa46-7bfb-4a1b-9fac-4ffb69e45277", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from scipy.special import erf\n", "import matplotlib.pyplot as plt\n", "import pyfvtool as pf" ] }, { "cell_type": "code", "execution_count": 2, "id": "ad3537ee-83db-4f81-b7f2-48638d6fd8bb", "metadata": {}, "outputs": [], "source": [ "# Spoiler alert!\n", "\n", "def T_analytic_dirichlet(x,t, alfa, T0, Ts):\n", " return (T0-Ts)*erf(x/np.sqrt(4*alfa*t))+Ts\n", "\n", "def T_analytic_neumann(x,t, alfa, T0, k, qs):\n", " return T0 + qs/k*np.sqrt(4*alfa*t/np.pi)*np.exp(-x**2/(4*alfa*t))\\\n", " - qs/k*x*(1-erf(x/np.sqrt(4*alfa*t)))" ] }, { "cell_type": "code", "execution_count": 3, "id": "df5f937e-551d-47b7-8604-0739f83b5758", "metadata": {}, "outputs": [], "source": [ "# Parameters\n", "L = 1.0 # [m] domain length\n", "k = 20.0 # 0.6 for water, 0.025 for air W/m/K\n", "rho = 8000.0 # kg/m^3\n", "c = 500.0 # J/kg/K (4200 for water, 1000 for air)\n", "alfa = k/(rho*c) # heat diffusion\n", "T0 = 300.0 # [K]\n", "Ts = 350.0 # [K]\n", "qs = 1000 # [W/m^2]\n", "t_sim = L**2/(20*alfa) # [s]\n", "time_steps = 50\n", "dt = t_sim/time_steps # \n", "Nx = 50 # number of cells" ] }, { "cell_type": "code", "execution_count": 4, "id": "ae1d51bc-4e13-43e9-88e9-3a7479f21610", "metadata": {}, "outputs": [], "source": [ "m = pf.Grid1D(Nx, L)" ] }, { "cell_type": "code", "execution_count": 5, "id": "f8f33dcd-ac71-4526-90d8-46eeeacde51e", "metadata": {}, "outputs": [], "source": [ "# Cell variable with initial condition\n", "T = pf.CellVariable(m, T0)" ] }, { "cell_type": "code", "execution_count": 6, "id": "cc6c5f88-467b-4f95-bcce-2fdbee91da4f", "metadata": {}, "outputs": [], "source": [ "# Boundary condition\n", "left_bc = \"Dirichlet\"\n", "\n", "BC = pf.BoundaryConditions(m)\n", "if left_bc == \"Dirichlet\":\n", " T.BCs.left.a[:] = 0.0\n", " T.BCs.left.b[:] = 1.0\n", " T.BCs.left.c[:] = Ts\n", " T_analytic = lambda x,t: T_analytic_dirichlet(x, t, alfa, T0, Ts)\n", "else:\n", " T.BCs.left.a[:] = k\n", " T.BCs.left.b[:] = 0.0\n", " T.BCs.left.c[:] = -qs\n", " T_analytic = lambda x,t: T_analytic_neumann(x, t, alfa, T0, k, qs)" ] }, { "cell_type": "code", "execution_count": 7, "id": "3d2335bc-6b3a-4098-a3b5-2196af34c23c", "metadata": {}, "outputs": [], "source": [ "# physical parameters\n", "alfa_cell = pf.CellVariable(m, alfa)\n", "alfa_face = pf.harmonicMean(alfa_cell)" ] }, { "cell_type": "code", "execution_count": 8, "id": "67534101-60dc-42f4-8624-2eeb40b1cc6d", "metadata": {}, "outputs": [], "source": [ "t=0\n", "while t" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.figure(1)\n", "plt.clf()\n", "plt.title('Left BC: '+left_bc)\n", "plt.plot(x, T_an, x, T_num, 'o')\n", "plt.legend({'Analytic', 'PyFVTool FVM'})\n", "plt.xlabel('x [m]')\n", "plt.ylabel('T [K]')\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 11, "id": "a56c0184-6fe6-4ae2-810b-a8fbe39f6e78", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.00021832771380885904\n" ] } ], "source": [ "er = np.sum(np.abs(T_num-T_an)/T_an)/Nx\n", "print(er)" ] }, { "cell_type": "code", "execution_count": null, "id": "7851ca90-a56e-40df-b536-b423455b6d45", "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.12.9" } }, "nbformat": 4, "nbformat_minor": 5 }