{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# SciPy" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "SciPy представляет собой совокупность математических алгоритмов и функций, построенных как расширение Numpy на Python. Он значительно расширяет возможности интерактивной сессии Python, предоставляя пользователю команды высокого уровня и классы для управления и визуализации данных. Со SciPy интерактивный сеанс Python становится средой обработки данных и системой прототипирования соперничающей с такими системами, как MATLAB, IDL, Octave, R-Lab, and SciLab." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*Модули SciPy*:\n", "\n", "* constants - Физические константы и коэффициенты пересчёта \n", "\n", "* cluster - Векторное квантование\n", "\n", "* fftpack - Дискретные алгоритмы преобразования Фурье \n", "\n", "* integrate - Инструменты для интегрирования \n", "\n", "* interpolate - Инструменты для интерполяции \n", "\n", "* io - Ввод-вывод данных \n", "\n", "* lib - Работа со сторонними библиотеками \n", "\n", "* linalg - Линейная алгебра \n", "\n", "* optimize - Средства оптимизации \n", "\n", "* sandbox - Экспериментальный код \n", "\n", "* signal - Обработка сигналов \n", "\n", "* sparse - Поддержка разреженных матриц \n", "\n", "* special - Специальные функции\n", "\n", "* stats - Статистические функции \n", "\n", "* weave - Использование кода, написанного на языках C и C++" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# подключим необходимые в дальнейшем модули\n", "import numpy as np\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Подмодуль linalg" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Модуль linalg предназначен для реализации линейной алгебры." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "from scipy import linalg " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Метод dot(a,b)**\n", "\n", "Матричное умножение a × b.\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A:\n", "[[1. 1. 1. 1.]\n", " [1. 1. 1. 1.]\n", " [1. 1. 1. 1.]\n", " [1. 1. 1. 1.]]\n", "B:\n", "[[1. 1.]\n", " [1. 1.]\n", " [1. 1.]\n", " [1. 1.]]\n", "[[4. 4.]\n", " [4. 4.]\n", " [4. 4.]\n", " [4. 4.]]\n" ] } ], "source": [ "A = np.array([[ 1., 1., 1., 1.], [ 1., 1., 1., 1.], [ 1., 1., 1., 1.],[ 1., 1., 1., 1.]])\n", "print('A:')\n", "print(A)\n", "B = np.array([[ 1., 1.], [ 1., 1.], [ 1., 1.], [ 1., 1.]])\n", "print('B:')\n", "print(B)\n", "\n", "q = np.dot(A, B)\n", "print(q)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Метод solve(a,b)**\n", "\n", "Задача: нужно найти такой x, что, если матрицу a умножить на x, то получится матрица b." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "a:\n", "[[ 3 2 0]\n", " [ 1 -1 0]\n", " [ 0 5 1]]\n", "b:\n", "[ 2 4 -1]\n", "x = [ 2. -2. 9.]\n" ] } ], "source": [ "a = np.array([[3, 2, 0],[1, -1, 0],[0, 5, 1]])\n", "print('a:')\n", "print(a)\n", "b = np.array([2, 4, -1])\n", "print('b:')\n", "print(b)\n", "x = linalg.solve(a,b)\n", "print('x =',x)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([ 2., 4., -1.])" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "np.dot(a, x) # Проверим - результат должен быть равен b" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Метод det(m)**\n", "\n", "Нахождение определителя матрицы." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[-4 -1 2]\n", " [10 4 -1]\n", " [ 8 3 1]]\n", "-14.000000000000009\n" ] } ], "source": [ "B = np.array([[-4, -1, 2], [10, 4, -1], [8, 3, 1]])\n", "print(B)\n", "q = np.linalg.det(B)\n", "print(q)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Метод svd(x)**\n", "\n", "Сингулярное разложение матрицы." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "x:\n", "[[-0.1472416 -0.24669754 0.61122716]\n", " [ 0.75507772 1.37491485 -0.50868729]\n", " [ 0.19862968 2.06707627 -0.25011693]\n", " [ 0.49005284 1.34922733 -0.77283106]]\n", "(4, 4) (3,) (3, 3)\n", " \n" ] } ], "source": [ "x = np.random.randn(4, 3)\n", "print('x:')\n", "print(x)\n", "U, D, V = linalg.svd(x)\n", "print (U.shape, D.shape, V.shape)\n", "print (type(U), type (D), type(V))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Подмодуль optimize" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "from scipy import optimize # подмодуль optimize" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "В модуле optimize есть **метод minimize**, которому достаточно передать функцию, передать начальное приближение. Мы видим, что значение функции равно примерно 2 в точке с координатами x: array([-1.03529084e-07, 4.99997747e+00]), которые с достаточно высокой степенью точности совпадают с правильными. " ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " fun: 2.000000274650647\n", " hess_inv: array([[ 4.54979753e-01, -8.89864272e-06],\n", " [-8.89864272e-06, 1.00000000e+00]])\n", " jac: array([-2.38418579e-07, 8.94069672e-08])\n", " message: 'Optimization terminated successfully.'\n", " nfev: 164\n", " nit: 36\n", " njev: 41\n", " status: 0\n", " success: True\n", " x: array([-1.03529084e-07, 4.99997747e+00])\n" ] } ], "source": [ "def f(x):\n", " return 3**(x[0]**2) + 3**(1e-8*x[1]**2)\n", "\n", "x_min = optimize.minimize(f, [-5, 5])\n", "print(x_min)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Если нужна только точка минимума, можно обратится к полю X нужного объекта" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[-1.03529084e-07 4.99997747e+00]\n" ] } ], "source": [ "print(x_min.x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Подмодуль integrate" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "import scipy.integrate as integrate # подмодуль integrate" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Функция quad** предназначена для интегрирования функции одной переменной между двумя точками. \n", "Точки могут быть (inf), чтобы указать бесконечные пределы. Например, предположим, что вы хотите интегрировать функцию \n", "Бесселя jv (2.5, x) вдоль интервала [0,4.5]" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(1.1178179380783253, 7.866317250224184e-09)\n" ] } ], "source": [ "import scipy.special as special\n", "result = integrate.quad(lambda x: special.jv(2.5,x), 0, 4.5)\n", "print(result)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Метод root**\n", "\n", "В следующем примере рассматривается трансцендентное уравнение с одной переменной\n", "x * x + 2cos (x) = 0.\n", "\n", "Корень уравнения можно найти следующим образом." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " fjac: array([[-1.]])\n", " fun: array([0.])\n", " message: 'The solution converged.'\n", " nfev: 10\n", " qtf: array([-2.77666778e-12])\n", " r: array([-3.3472241])\n", " status: 1\n", " success: True\n", " x: array([-0.73908513])\n" ] } ], "source": [ "from scipy.optimize import root\n", "def func(x):\n", " return x*2 + 2 * np.cos(x)\n", "sol = root(func, 0.3)\n", "print(sol)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## sparse - Поддержка разреженных матриц" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Работа с большими разреженными матрицами в Python: десятки миллионов элементов = 10^4 x 10^4. Пример - исходная матрица не помещается в ОЗУ:\n", "\n", "https://docs.scipy.org/doc/scipy/reference/sparse.html\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from scipy.sparse import lil_matrix\n", "from scipy.sparse.linalg import spsolve\n", "from numpy.linalg import solve, norm\n", "from numpy.random import rand" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Создадим матрицу 1000x1000 lil_matrix (lil = list of lists):" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "A = lil_matrix((1000, 1000))\n", "A[0, :100] = rand(100)\n", "A[1, 100:200] = A[0, :100]\n", "A.setdiag(rand(1000))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Преобразуем ее к формату CSR и решим A*x = b:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "A = A.tocsr()\n", "b = rand(1000)\n", "x = spsolve(A, b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Преобразуем матрицу к неразреженной, решим, и проверим, что результаты совпадают:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "x_ = solve(A.toarray(), b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Теперь вычислим норму ошибки (ответы почти не отличаются):" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "err = norm(x-x_)\n", "err < 1e-10" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "from scipy.sparse import rand\n", "import numpy\n", "from datetime import datetime\n", "\n", "n = 10000\n", "A = rand(n, n, density=0.25, format=\"csr\")\n", "b = numpy.random.rand(n)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "48" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import sys\n", " \n", "sys.getsizeof(A)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "80112" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sys.getsizeof(b)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0:11:59.469898\n" ] } ], "source": [ "start_time = datetime.now()\n", "x = spsolve(A, b)\n", "print(datetime.now() - start_time)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "A = A.toarray()" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "800000128" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sys.getsizeof(A)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0:00:38.450504\n" ] } ], "source": [ "start_time = datetime.now()\n", "x_ = solve(A, b)\n", "print(datetime.now() - start_time)" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "err = norm(x-x_)\n", "err < 1e-8" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "32" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sys.getsizeof(1111111111111111134)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755626288" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x=3\n", "id(x)" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755626288" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "y=3\n", "id(y)" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755626352" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x=5\n", "id(x)" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755626352" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "v=1+4\n", "id(v)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755634384" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x=256\n", "id(x)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1512755634384" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "y=256\n", "id(y)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1512755626032\n", "1512755626032\n" ] } ], "source": [ "x=-5\n", "print(id(x))\n", "y=-5\n", "print(id(y))" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1512808410480\n", "1512808411632\n" ] } ], "source": [ "x=257\n", "print(id(x))\n", "y=257\n", "print(id(y))" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1512808410512\n", "1512808415120\n" ] } ], "source": [ "x=-6\n", "print(id(x))\n", "y=-6\n", "print(id(y))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "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.10.9" } }, "nbformat": 4, "nbformat_minor": 4 }