C05 - count (recurrence).ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Первые двадцать = 0 1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597 2584 4181\n",
      "Через 5 часов = 514229\n"
     ]
    }
   ],
   "source": [
    "import math\n",
    "import numpy as np\n",
    "total_time = 5 * 60 // 10  # число 10-минуток в 5 часах\n",
    "\n",
    "# для начала приведём решение для выражения n[i] = n[i-1] + n[i-2]\n",
    "values = [0, 1] # начальные условия для вычисления рекуррентного выражения\n",
    "\n",
    "# сгеренируем эти числа явным образом и проверим:\n",
    "def next1(prev2, prev1):\n",
    "    return prev1 + prev2\n",
    "\n",
    "for i in range(total_time - 2):\n",
    "    values.append(next1(*values[-2:]))\n",
    "\n",
    "print(\"Первые двадцать =\", *values[:20])\n",
    "print(\"Через 5 часов =\", values[-1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "r1, r2 = -0.618, 1.618\n"
     ]
    }
   ],
   "source": [
    "# теперь построим функцию для нахождения членов последовательности в явном виде\n",
    "# характеристическое уравнение\n",
    "# r**2 - r - 1 = 0\n",
    "a, b, c = 1, -1, -1\n",
    "# решаем через дискриминант, как в школе\n",
    "D = b ** 2 - 4 * a * c\n",
    "r1, r2 = (-b - math.sqrt(D)) / (2 * a), (-b + math.sqrt(D)) / (2 * a)\n",
    "print(\"r1, r2 = {:.3f}, {:.3f}\".format(r1, r2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "α, β = -0.447, 0.447\n",
      "Первые двадцать = 0 1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597 2584 4181\n",
      "Через 5 часов = 514229\n",
      "Сейчас заражены и активны 832040 машин\n"
     ]
    }
   ],
   "source": [
    "# теперь решим систему относительно α и β\n",
    "# α*r1^0 + β*r1^0 = n[0]\n",
    "# α*r1^1 + β*r2^1 = n[1]\n",
    "\n",
    "# Решаем явно через обращение матрицы.  A*x=y -> x=(A^-1)*y\n",
    "A = np.matrix([[r1**i, r2**i] for i in range(2)])\n",
    "y = np.matrix([values[:2]]).T\n",
    "\n",
    "x = A.I * y\n",
    "α, β = x.T.A[0]\n",
    "print(\"α, β = {:.3f}, {:.3f}\".format(α, β))\n",
    "\n",
    "# проверяем построенное решение на корректность\n",
    "def gen1(i):\n",
    "    return int(round(α * (r1 ** i) + β * (r2 ** i)))\n",
    "\n",
    "print(\"Первые двадцать =\", *[gen1(i) for i in range(20)])\n",
    "print(\"Через 5 часов =\", gen1(total_time - 1))\n",
    "\n",
    "# итого, число активных машин к моменту времени t равно\n",
    "t = total_time\n",
    "result = gen1(t - 1) + gen1(t - 2)\n",
    "print(\"Сейчас заражены и активны {} машин\".format(result))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Первые двадцать = 0 1 2 4 7 12 20 33 54 88 143 232 376 609 986 1596 2583 4180 6764 10945\n",
      "Через 5 часов паролей украдено = 1346268\n"
     ]
    }
   ],
   "source": [
    "# теперь приведём решение для выражения n[i] = 2*n[i-1] - n[i-3]\n",
    "def next2(prev3, prev2, prev1):\n",
    "    return 2 * prev1 - prev3\n",
    "\n",
    "compr_values = [0, 1, 2]   # начальные значения\n",
    "\n",
    "# сгеренируем эти числа явным образом\n",
    "for i in range(total_time - 3):\n",
    "    compr_values.append(next2(*compr_values[-3:]))\n",
    "\n",
    "print(\"Первые двадцать =\", *compr_values[:20])\n",
    "print(\"Через 5 часов паролей украдено =\", compr_values[-1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "r1, r2, r3 = -0.618, 1.000, 1.618\n",
      "α, β, γ == -0.171, -1.000, 1.171\n",
      "Первые двадцать = 0 1 2 4 7 12 20 33 54 88 143 232 376 609 986 1596 2583 4180 6764 10945\n",
      "Через 5 часов паролей украдено = 1346268\n"
     ]
    }
   ],
   "source": [
    "# теперь построим решение в явном виде\n",
    "\n",
    "# характеристическое уравнение\n",
    "# r^3 - 2*r^2 - 0*r + 1 = 0\n",
    "# (r-1) * (r^2 - r - 1) = 0 - решим вторую скобку\n",
    "a, b, c = 1, -1, -1\n",
    "D = b ** 2 - 4 * a * c\n",
    "r1, r2, r3 = (-b - math.sqrt(D)) / (2 * a), 1.0, (-b + math.sqrt(D)) / (2 * a)\n",
    "print(\"r1, r2, r3 = {:.3f}, {:.3f}, {:.3f}\".format(r1, r2, r3))\n",
    "\n",
    "# найдём коэффициенты\n",
    "A = np.matrix([[r1**i, r2**i, r3**i] for i in range(3)])\n",
    "y = np.matrix([compr_values[:3]]).T\n",
    "r = A.I * y\n",
    "\n",
    "α, β, γ = r.T.A[0]\n",
    "print(\"α, β, γ == {:.3f}, {:.3f}, {:.3f}\".format(α, β, γ))\n",
    "\n",
    "def gen2(i):\n",
    "    return int(round(α*(r1**i)+β*(r2**i)+γ*(r3**i)))\n",
    "\n",
    "print(\"Первые двадцать =\", *[gen2(i) for i in range(20)])\n",
    "print(\"Через 5 часов паролей украдено =\", gen2(total_time - 1))"
   ]
  }
 ],
 "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.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}