{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "fda98006",
   "metadata": {},
   "source": [
    "# Solution of second-order linear ordinary differential equation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a4b21c72",
   "metadata": {},
   "source": [
    "## Input vars"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6653253d",
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 20\n",
    "linear_file = f\"matrix_linear_{N}.txt\"\n",
    "cubic_file = f\"matrix_cubic_{N}.txt\""
   ]
  },
  {
   "cell_type": "markdown",
   "id": "781047cc",
   "metadata": {},
   "source": [
    "## System cells"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7fced95c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def write_to_file(name: str, array):\n",
    "    with open(f\"{name}\", 'w') as file:\n",
    "        for element in array:\n",
    "            file.write(f\"{element} \")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d20bde93",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8af07a56",
   "metadata": {},
   "source": [
    "### Solution of $5u'' + 4u' + 1 = 0, u'(0) = u(0), u(10) = 5$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d0fe0703",
   "metadata": {},
   "outputs": [],
   "source": [
    "def u(x):\n",
    "    return (9 * np.exp(8) * (x - 30) + 155 * np.exp(8 - 4 * x / 5) - 5 * (x + 1)) / (20 - 36 * np.exp(8))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "efac514a",
   "metadata": {},
   "source": [
    "### Linear element"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "18bd929b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def show_plot_linear(filename: str, elements: int):\n",
    "    with open(filename, 'r') as file:\n",
    "        data_str = file.read()\n",
    "\n",
    "    x = np.linspace(0, 10, elements + 1)\n",
    "    y = np.fromstring(data_str, sep=' ')\n",
    "    y_real = u(x)\n",
    "\n",
    "    plt.plot(x, y, label=\"Linear element solution\", color='blue')\n",
    "    plt.plot(x, y_real, label=\"Numeral solution\", color='black')\n",
    "    plt.title(f\"Linear element, elements = {elements}\")\n",
    "    plt.grid(True)\n",
    "    plt.xlabel(\"x\")\n",
    "    plt.ylabel(\"u(x)\")\n",
    "    plt.legend()\n",
    "    plt.savefig(f\"linear_{elements}.png\", dpi=300)\n",
    "    plt.show()\n",
    "\n",
    "    write_to_file(f\"linear_real_y_{elements}.txt\", y_real)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d9a38740",
   "metadata": {},
   "source": [
    "### Cubic element"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "73525bc7",
   "metadata": {},
   "outputs": [],
   "source": [
    "def show_plot_cubic(filename: str, elements: int):\n",
    "    with open(filename, 'r') as file:\n",
    "        data_str = file.read()\n",
    "\n",
    "    x = np.linspace(0, 10, 3 * elements + 1)\n",
    "    y = np.fromstring(data_str, sep=' ')\n",
    "    y_real = u(x)\n",
    "\n",
    "    plt.plot(x, y, label=\"Cubic element solution\", color=\"red\")\n",
    "    plt.plot(x, y_real, label=\"Numeral solution\", color='black')\n",
    "    plt.title(f\"Cubic element, elements = {elements} \")\n",
    "    plt.xlabel(\"x\")\n",
    "    plt.ylabel(\"u(x)\")\n",
    "    plt.grid(True)\n",
    "    plt.legend()\n",
    "    plt.savefig(f\"cubic_{elements}.png\", dpi=300)\n",
    "    plt.show()\n",
    "\n",
    "    write_to_file(f\"cubic_real_y_{elements}.txt\", y_real)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3584baa4",
   "metadata": {},
   "source": [
    "## Main part"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b8a772e3",
   "metadata": {},
   "outputs": [],
   "source": [
    "show_plot_linear(linear_file, N)\n",
    "show_plot_cubic(cubic_file, N)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22fe3de6",
   "metadata": {},
   "source": [
    "## Errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "31de1104",
   "metadata": {},
   "outputs": [],
   "source": [
    "def read_data_from_file(filename):\n",
    "    with open(filename, 'r') as file:\n",
    "        content = file.read().strip()\n",
    "        data = np.array([float(x) for x in content.split()])\n",
    "    return data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e0e9266d",
   "metadata": {},
   "outputs": [],
   "source": [
    "def show_error(program_file: str, analytical_file: str):\n",
    "    # Считывание данных из файлов\n",
    "    real = read_data_from_file(analytical_file)\n",
    "    matrix = read_data_from_file(program_file)\n",
    "    \n",
    "    absolute_errors = np.abs(real - matrix)\n",
    "    relative_errors = absolute_errors / (np.abs(real))\n",
    "    relative_errors_percent = relative_errors * 100\n",
    "\n",
    "    max_relative_error = np.max(relative_errors)\n",
    "    max_relative_error_percent = np.max(relative_errors_percent)\n",
    "    max_error_index = np.argmax(relative_errors)\n",
    "\n",
    "    print(f\"\\nМаксимальная относительная погрешность: {max_relative_error:.6f} ({max_relative_error_percent:.2f}%)\")\n",
    "    print(f\"Индекс точки с максимальной погрешностью: {max_error_index}\")\n",
    "    print(f\"Значения в точке с максимальной погрешностью:\")\n",
    "    print(f\"  cubic_real: {real[max_error_index]}\")\n",
    "    print(f\"  matrix_cubic: {matrix[max_error_index]}\")\n",
    "    print(f\"  Абсолютная разность: {absolute_errors[max_error_index]:.6f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "89918171",
   "metadata": {},
   "outputs": [],
   "source": [
    "show_error(linear_file, f\"linear_real_y_{N}.txt\")\n",
    "show_error(cubic_file, f\"cubic_real_y_{N}.txt\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}