{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "code",
      "source": [
        "import math\n",
        "import time\n",
        "from pathlib import Path"
      ],
      "metadata": {
        "id": "C0KZdYS6MZvR"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "L = 0.02 #размер трубы в м\n",
        "num_cell = 100 #число ячеек\n",
        "\n",
        "dx = L / num_cell\n",
        "\n",
        "n_2 = 4e11 #сколько реальных частиц моделирует одна одномерная частица\n",
        "\n",
        "N_part = num_cell * 50\n",
        "\n",
        "T = 3.5e-9 #время расчета\n",
        "dt = 1e-11 #шаг по времени\n",
        "\n",
        "eps0 = 8.85e-12\n",
        "\n",
        "denom = dx * dx * eps0"
      ],
      "metadata": {
        "id": "53AVcwzi49ha"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "E = [0 for _ in range(num_cell + 1)]\n",
        "fi = [0 for _ in range(num_cell + 1)]\n",
        "\n",
        "p = [0 for _ in range(num_cell - 1)]\n",
        "q = [0 for _ in range(num_cell - 1)]"
      ],
      "metadata": {
        "id": "xeiRwSLmCbi-"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "def solve_field(rho):\n",
        "  fi[0] = 0\n",
        "  fi[-1] = 0\n",
        "\n",
        "  p[1] = 0.5\n",
        "  q[1] = (dx * dx * rho[1] - fi[0]) / (-2)\n",
        "\n",
        "  for i in range(2, len(fi) - 2, 1):\n",
        "    p[i] = -1 / (p[i - 1] - 2)\n",
        "    q[i] = (dx * dx * rho[i] - q[i-1]) / (p[i-1] - 2)\n",
        "  fi[-2] = (dx * dx * rho[-2] - fi[-1] - q[-3]) / (p[-3] - 2)\n",
        "\n",
        "  for i in range(len(fi) - 3, 0, -1):\n",
        "    fi[i] = p[i] * fi[i + 1] + q[i]\n",
        "\n",
        "  E[0] = (fi[0] - fi[1]) / dx\n",
        "  for i in range(1, num_cell - 1):\n",
        "    E[i] = (fi[i-1] - fi[i+1]) / (2 * dx)\n",
        "  E[-1] = (fi[-2] - fi[-1]) / dx"
      ],
      "metadata": {
        "id": "ggskSz2KCRxc"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "YkhrWNtDr9cn"
      },
      "outputs": [],
      "source": [
        " class Particles:\n",
        "  def __init__(self, k=4, ioe=\"+\", T_ev = 10, num_cell=num_cell):\n",
        "    self.q = (1.6e-19 if ioe == '+' else -1.6e-19)\n",
        "    self.m = (1.67e-27 if ioe == '+' else 9.31e-31)\n",
        "    self.v = [0 for _ in range(N_part)]\n",
        "    self.rho = [0 for _ in range(num_cell + 1)]\n",
        "\n",
        "    self.x = [i / N_part * L for i in range(N_part)]\n",
        "\n",
        "    if ioe == '-':\n",
        "      for i in range(N_part):\n",
        "\n",
        "        self.x[i] += 0.0005 * L * math.sin(self.x[i] / L * 2 * math.pi * k)\n",
        "\n",
        "  def move(self, E):\n",
        "    for i in range(N_part):\n",
        "      nc = math.floor(self.x[i] // dx)\n",
        "      x_loc = self.x[i] % dx\n",
        "      E_by_x = (E[nc] * (dx - x_loc) + E[nc + 1] * x_loc) / dx\n",
        "      self.v[i] += E_by_x * self.q / self.m * dt\n",
        "\n",
        "      self.x[i] += self.v[i] * dt\n",
        "\n",
        "      #ГУ при выходе за границу\n",
        "      if self.x[i] > L:\n",
        "        self.x[i] -= 2 * (self.x[i] - L)\n",
        "      elif self.x[i] == L:\n",
        "        self.x[i] *= 0.9995\n",
        "      elif self.x[i] < 0:\n",
        "        self.x[i] *= (-1)\n",
        "\n",
        "  def CIC(self):\n",
        "    for i in range(len(self.rho)):\n",
        "      self.rho[i] = 0\n",
        "\n",
        "    for i in range(len(self.x)):\n",
        "\n",
        "      nc = int(self.x[i] // dx)\n",
        "      #print(nc)\n",
        "      x_loc = self.x[i] % dx\n",
        "\n",
        "      self.rho[nc] += self.q * (dx - x_loc) * (-n_2) / denom\n",
        "      self.rho[nc + 1] += self.q * x_loc * (-n_2) / denom\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "el = Particles(ioe='-')\n",
        "ions = Particles()\n"
      ],
      "metadata": {
        "id": "7n6jF_ym8P8x"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "!pwd"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "LrzEVqrjBTT7",
        "outputId": "67aef206-24e0-42be-caa0-d5e62bbbac82"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "/content\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "start = time.time()\n",
        "\n",
        "for j in range(1, int(T // dt) + 1):\n",
        "  el.CIC()\n",
        "  ions.CIC()\n",
        "  #print('плотность расчиталась')\n",
        "\n",
        "  rho = [el.rho[i] + ions.rho[i] for i in range(num_cell)]\n",
        "  solve_field(rho)\n",
        "\n",
        "\n",
        "  el.move(E)\n",
        "  ions.move(E)\n",
        "\n",
        "  if True: #j % 1000 == 0:\n",
        "\n",
        "    rho_path = Path(\"rho.txt\")\n",
        "\n",
        "    with rho_path.open(\"a\", encoding=\"utf-8\") as f:\n",
        "      for item in el.rho:\n",
        "        f.write(f\"{item} \")\n",
        "      f.write('\\n')\n",
        "\n",
        "end = time.time()\n",
        "\n",
        "print(f\"Время выполнения: {end - start} секунд\")\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "xSTjBa_IN7V-",
        "outputId": "2f5680ce-dbb0-48d3-c5e1-802cc73d222e"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Время выполнения: 3.641313314437866 секунд\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "import math\n",
        "import cv2\n",
        "import os"
      ],
      "metadata": {
        "id": "y61eZDBgCuNk"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "x_ions = []\n",
        "with open ('rho.txt') as data:\n",
        "    for line in data:\n",
        "        x_ions.append(line.split())\n",
        "        #x_ions = x_ions.astype(np.float64)\n",
        "    #print(type(data))\n",
        "x_ions = np.array(x_ions).astype(np.float64)"
      ],
      "metadata": {
        "id": "J-9-W6oKC0ec"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "!rm rho_*.png"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Cj0XjTR6FL2M",
        "outputId": "0cf739ea-4acc-4459-b625-306fe65c09e5"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "rm: cannot remove 'rho_*.png': No such file or directory\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "i = 0\n",
        "for i in range(len(x_ions)):\n",
        "    plt.figure()\n",
        "    plt.title('rho')\n",
        "    plt.ylim(x_ions.min(), x_ions.max())\n",
        "    #plt.bar(np.arange(len(x_el[0])), x_el[i], alpha = 0.3)\n",
        "    plt.bar(np.arange(len(x_ions[0])), x_ions[i])\n",
        "    plt.savefig(f'rho_{i}')\n",
        "    plt.close()"
      ],
      "metadata": {
        "id": "EBPUD4TrC6Ib"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "image_folder = '.'  # Папка с изображениями (текущая папка)\n",
        "output_video = '!_rho.avi'  # Имя выходного видеофайла\n",
        "fps = 5  # Частота кадров\n",
        "\n",
        "image_files = [f for f in os.listdir(image_folder) if f.startswith('rho_') and f.endswith('.png')]\n",
        "\n",
        "image_files.sort(key=lambda f: int(f.split('_')[1].split('.')[0]))\n",
        "\n",
        "first_image = cv2.imread(os.path.join(image_folder, image_files[0]))\n",
        "height, width, _ = first_image.shape\n",
        "\n",
        "fourcc = cv2.VideoWriter_fourcc(*'XVID')\n",
        "video_writer = cv2.VideoWriter(output_video, fourcc, fps, (width, height))\n",
        "\n",
        "for image_file in image_files:\n",
        "    image_path = os.path.join(image_folder, image_file)\n",
        "    image = cv2.imread(image_path)\n",
        "    video_writer.write(image)\n",
        "\n",
        "video_writer.release()\n",
        "\n",
        "print(f\"Видео успешно создано: {output_video}\")"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "d5ao6mujDu9l",
        "outputId": "aaba809c-e351-41ea-9402-e94b0845fcdd"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Видео успешно создано: !_rho.avi\n"
          ]
        }
      ]
    }
  ]
}