Last active
November 3, 2016 22:49
-
-
Save nmayorov/a16f3d9e90cfc7c86a96b1bc9e68b19f to your computer and use it in GitHub Desktop.
Profiling of ODE solvers
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
{ | |
"cells": [ | |
{ | |
"cell_type": "code", | |
"execution_count": 1, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [], | |
"source": [ | |
"%load_ext line_profiler" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 2, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"import numpy as np" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 3, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"from scipy.integrate import solve_ivp" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 4, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"def fun(t, u, m):\n", | |
" N = u.shape[0] // 4\n", | |
" \n", | |
" x = u[:N]\n", | |
" y = u[N:2*N]\n", | |
" x_dot = u[2*N:3*N]\n", | |
" y_dot = u[3*N:]\n", | |
" \n", | |
" dx = x[:, None] - x\n", | |
" dy = y[:, None] - y\n", | |
" denom = (dx ** 2 + dy ** 2) ** 1.5\n", | |
" i = np.arange(N)\n", | |
" denom[i, i] = 1\n", | |
"\n", | |
" Fx = m * dx / denom\n", | |
" Fx[i, i] = 0\n", | |
" Fx = -np.sum(Fx, axis=1)\n", | |
" \n", | |
" Fy = m * dy / denom\n", | |
" Fy[i, i] = 0\n", | |
" Fy = -np.sum(Fy, axis=1)\n", | |
" \n", | |
" return np.hstack((x_dot, y_dot, Fx, Fy))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 5, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"def jac(t, u, m):\n", | |
" N = u.shape[0] // 4\n", | |
" \n", | |
" x = u[:N]\n", | |
" y = u[N:2*N]\n", | |
" x_dot = u[2*N:3*N]\n", | |
" y_dot = u[3*N:]\n", | |
" \n", | |
" J = np.zeros((4 * N, 4 * N))\n", | |
" i = np.arange(N)\n", | |
" J[i, 2 * N + i] = 1\n", | |
" J[i + N, 3 * N + i] = 1\n", | |
" \n", | |
" dx = x[:, None] - x\n", | |
" dy = y[:, None] - y\n", | |
"\n", | |
" r = (dx ** 2 + dy ** 2) ** 0.5\n", | |
" r[i, i] = 1\n", | |
" r_m3 = r ** -3\n", | |
" r_m3[i, i] = 0\n", | |
" r_m5 = r ** -5\n", | |
" r_m5[i, i] = 0\n", | |
" \n", | |
" J[2*N:3*N, :N] = m * (r_m3 - 3 * dx**2 * r_m5)\n", | |
" J[2*N:3*N, N:2*N] = -3 * m * dx * dy * r_m5\n", | |
"\n", | |
" J[2*N + i, i] = -np.dot(r_m3, m) + 3 * np.dot(dx**2 * r_m5, m)\n", | |
" J[2*N + i, N + i] = 3 * np.dot(dx * dy * r_m5, m)\n", | |
" \n", | |
" J[3*N:, :N] = -3 * m * dx * dy * r_m5\n", | |
" J[3*N:, N:2*N] = m * (r_m3 - 3 * dy**2 * r_m5)\n", | |
" J[3*N + i, i] = 3 * np.dot(dx * dy * r_m5, m)\n", | |
" J[3*N + i, N + i] = -np.dot(r_m3, m) + 3 * np.dot(dy**2 * r_m5, m)\n", | |
" \n", | |
" \n", | |
" return J" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 6, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"def setup(N, m0=1000, r0=1):\n", | |
" x0 = np.zeros(N)\n", | |
" y0 = np.zeros(N)\n", | |
" x_dot_0 = np.zeros(N)\n", | |
" y_dot_0 = np.zeros(N)\n", | |
" m = np.zeros(N)\n", | |
" \n", | |
" m[0] = m0\n", | |
" m[1:] = 1 \n", | |
" for i in range(1, N):\n", | |
" r = i * r0\n", | |
" \n", | |
" phi = np.random.uniform(-np.pi, np.pi)\n", | |
" x0[i] = r * np.cos(phi)\n", | |
" y0[i] = r * np.sin(phi)\n", | |
"\n", | |
" omega = (-1) ** i * (m[0] / r ** 3) ** 0.5\n", | |
" \n", | |
" x_dot_0[i] = -r * omega * np.sin(phi)\n", | |
" y_dot_0[i] = r * omega * np.cos(phi)\n", | |
" \n", | |
" return np.hstack((x0, y0, x_dot_0, y_dot_0)), m " | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"# RK methods" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Small N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 8, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 3" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 9, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 20]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 10, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 11, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 12, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 593042 function calls in 1.058 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 50 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 12752 0.400 0.000 0.770 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 2125 0.099 0.000 0.945 0.000 rk.py:15(rk_step)\n", | |
" 25505 0.085 0.000 0.085 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 51008 0.071 0.000 0.138 0.000 shape_base.py:9(atleast_1d)\n", | |
" 17003 0.054 0.000 0.054 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 1596 0.041 0.000 1.017 0.001 rk.py:104(_step_impl)\n", | |
" 52605 0.038 0.000 0.047 0.000 numeric.py:484(asanyarray)\n", | |
" 25504 0.033 0.000 0.137 0.000 fromnumeric.py:1743(sum)\n", | |
" 12753 0.030 0.000 0.030 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 12752 0.027 0.000 0.164 0.000 shape_base.py:275(<listcomp>)\n", | |
" 12752 0.020 0.000 0.020 0.000 {built-in method numpy.core.multiarray.arange}\n", | |
" 12752 0.018 0.000 0.212 0.000 shape_base.py:232(hstack)\n", | |
" 12752 0.013 0.000 0.820 0.000 base.py:124(fun)\n", | |
" 12752 0.013 0.000 0.783 0.000 <string>:1(<lambda>)\n", | |
" 12752 0.013 0.000 0.807 0.000 base.py:11(fun_wrapped)\n", | |
" 2128 0.011 0.000 0.047 0.000 linalg.py:1976(norm)\n", | |
" 67488 0.011 0.000 0.011 0.000 {built-in method numpy.core.multiarray.array}\n", | |
" 14882 0.011 0.000 0.013 0.000 numeric.py:414(asarray)\n", | |
" 55794 0.011 0.000 0.011 0.000 {method 'append' of 'list' objects}\n", | |
" 106804 0.011 0.000 0.011 0.000 {built-in method builtins.len}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10bf9f710>" | |
] | |
}, | |
"execution_count": 12, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Medium N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 13, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 50" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 14, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 10]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 15, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 16, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 17, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 304256 function calls in 1.341 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 50 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 6578 0.922 0.000 1.165 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 13157 0.077 0.000 0.077 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 1096 0.069 0.000 1.287 0.001 rk.py:15(rk_step)\n", | |
" 26312 0.040 0.000 0.076 0.000 shape_base.py:9(atleast_1d)\n", | |
" 8771 0.026 0.000 0.026 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 683 0.026 0.000 1.331 0.002 rk.py:104(_step_impl)\n", | |
" 26996 0.020 0.000 0.025 0.000 numeric.py:484(asanyarray)\n", | |
" 13156 0.020 0.000 0.109 0.000 fromnumeric.py:1743(sum)\n", | |
" 6579 0.019 0.000 0.019 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 6578 0.015 0.000 0.091 0.000 shape_base.py:275(<listcomp>)\n", | |
" 6578 0.013 0.000 0.013 0.000 {built-in method numpy.core.multiarray.arange}\n", | |
" 6578 0.011 0.000 0.120 0.000 shape_base.py:232(hstack)\n", | |
" 6578 0.009 0.000 1.174 0.000 <string>:1(<lambda>)\n", | |
" 6578 0.008 0.000 1.195 0.000 base.py:124(fun)\n", | |
" 6578 0.007 0.000 1.187 0.000 base.py:11(fun_wrapped)\n", | |
" 13156 0.006 0.000 0.084 0.000 _methods.py:31(_sum)\n", | |
" 28359 0.006 0.000 0.006 0.000 {method 'append' of 'list' objects}\n", | |
" 1099 0.006 0.000 0.014 0.000 linalg.py:1976(norm)\n", | |
" 7679 0.006 0.000 0.007 0.000 numeric.py:414(asarray)\n", | |
" 34676 0.006 0.000 0.006 0.000 {built-in method numpy.core.multiarray.array}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10bfa7748>" | |
] | |
}, | |
"execution_count": 17, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": { | |
"collapsed": true | |
}, | |
"source": [ | |
"## Large N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 22, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 200" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 23, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 1]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 24, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 25, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 26, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 248576 function calls in 10.451 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 50 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 5378 9.483 0.002 10.091 0.002 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 10757 0.322 0.000 0.322 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 896 0.132 0.000 10.365 0.012 rk.py:15(rk_step)\n", | |
" 7171 0.085 0.000 0.085 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 21512 0.057 0.000 0.103 0.000 shape_base.py:9(atleast_1d)\n", | |
" 10756 0.048 0.000 0.395 0.000 fromnumeric.py:1743(sum)\n", | |
" 543 0.040 0.000 10.431 0.019 rk.py:104(_step_impl)\n", | |
" 5378 0.035 0.000 0.035 0.000 {built-in method numpy.core.multiarray.arange}\n", | |
" 5379 0.034 0.000 0.034 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 5378 0.029 0.000 10.120 0.002 <string>:1(<lambda>)\n", | |
" 22056 0.027 0.000 0.033 0.000 numeric.py:484(asanyarray)\n", | |
" 5378 0.023 0.000 0.126 0.000 shape_base.py:275(<listcomp>)\n", | |
" 5378 0.023 0.000 0.178 0.000 shape_base.py:232(hstack)\n", | |
" 5378 0.016 0.000 10.157 0.002 base.py:124(fun)\n", | |
" 10756 0.014 0.000 0.336 0.000 _methods.py:31(_sum)\n", | |
" 5378 0.013 0.000 10.140 0.002 base.py:11(fun_wrapped)\n", | |
" 10756 0.011 0.000 0.011 0.000 {built-in method builtins.isinstance}\n", | |
" 899 0.010 0.000 0.020 0.000 linalg.py:1976(norm)\n", | |
" 6279 0.008 0.000 0.009 0.000 numeric.py:414(asarray)\n", | |
" 28336 0.007 0.000 0.007 0.000 {built-in method numpy.core.multiarray.array}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10bf9ffd0>" | |
] | |
}, | |
"execution_count": 26, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"# Radau method" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Small N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 27, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 3" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 28, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 20]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 29, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 31, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='Radau', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 32, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 1762636 function calls in 3.066 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 102 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 17868 0.601 0.000 1.167 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 2334 0.269 0.000 2.203 0.001 radau.py:49(solve_collocation_system)\n", | |
" 33096 0.208 0.000 0.276 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 57723 0.208 0.000 0.208 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 1552 0.136 0.000 3.060 0.002 radau.py:381(_step_impl)\n", | |
" 16548 0.110 0.000 0.573 0.000 blas.py:226(_get_funcs)\n", | |
" 71472 0.109 0.000 0.207 0.000 shape_base.py:9(atleast_1d)\n", | |
" 1071 0.106 0.000 0.118 0.000 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 12438 0.091 0.000 0.704 0.000 decomp_lu.py:75(lu_solve)\n", | |
" 16548 0.076 0.000 0.428 0.000 blas.py:177(find_best_blas_type)\n", | |
" 16548 0.062 0.000 0.164 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 19438 0.060 0.000 0.060 0.000 {method 'dot' of 'numpy.ndarray' objects}\n", | |
" 78464 0.060 0.000 0.074 0.000 numeric.py:484(asanyarray)\n", | |
" 35736 0.053 0.000 0.210 0.000 fromnumeric.py:1743(sum)\n", | |
" 190446 0.049 0.000 0.049 0.000 numerictypes.py:951(<listcomp>)\n", | |
" 4110 0.048 0.000 0.187 0.000 decomp_lu.py:17(lu_factor)\n", | |
" 17869 0.047 0.000 0.047 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 7003 0.045 0.000 0.088 0.000 linalg.py:1976(norm)\n", | |
" 17868 0.042 0.000 0.248 0.000 shape_base.py:275(<listcomp>)\n", | |
" 44504 0.039 0.000 0.051 0.000 numeric.py:414(asarray)\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10bfa3630>" | |
] | |
}, | |
"execution_count": 32, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Medium N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 33, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 20" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 34, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 10]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 35, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 37, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='Radau', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 38, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 1019039 function calls in 2.586 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 102 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 10288 0.515 0.000 0.845 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 2318 0.405 0.000 0.551 0.000 decomp_lu.py:17(lu_factor)\n", | |
" 1436 0.170 0.000 1.561 0.001 radau.py:49(solve_collocation_system)\n", | |
" 7138 0.156 0.000 0.502 0.000 decomp_lu.py:75(lu_solve)\n", | |
" 818 0.143 0.000 2.576 0.003 radau.py:381(_step_impl)\n", | |
" 33189 0.135 0.000 0.135 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 18912 0.116 0.000 0.153 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 9456 0.095 0.000 0.158 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 516 0.068 0.000 0.075 0.000 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 9456 0.064 0.000 0.325 0.000 blas.py:226(_get_funcs)\n", | |
" 41152 0.061 0.000 0.115 0.000 shape_base.py:9(atleast_1d)\n", | |
" 11120 0.045 0.000 0.045 0.000 {method 'dot' of 'numpy.ndarray' objects}\n", | |
" 9456 0.043 0.000 0.240 0.000 blas.py:177(find_best_blas_type)\n", | |
" 45128 0.033 0.000 0.041 0.000 numeric.py:484(asanyarray)\n", | |
" 20576 0.030 0.000 0.130 0.000 fromnumeric.py:1743(sum)\n", | |
" 109400 0.027 0.000 0.027 0.000 numerictypes.py:951(<listcomp>)\n", | |
" 10289 0.027 0.000 0.027 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 3985 0.024 0.000 0.049 0.000 linalg.py:1976(norm)\n", | |
" 10288 0.023 0.000 0.138 0.000 shape_base.py:275(<listcomp>)\n", | |
" 9093 0.023 0.000 0.023 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10bfc4b38>" | |
] | |
}, | |
"execution_count": 38, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Large N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 48, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 100" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 49, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 1]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 50, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 51, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='Radau', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 53, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 91912 function calls in 3.728 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 102 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 218 1.951 0.009 2.137 0.010 decomp_lu.py:17(lu_factor)\n", | |
" 929 0.484 0.001 0.551 0.001 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 69 0.452 0.007 3.722 0.054 radau.py:381(_step_impl)\n", | |
" 641 0.339 0.001 0.398 0.001 decomp_lu.py:75(lu_solve)\n", | |
" 859 0.169 0.000 0.188 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 52 0.067 0.001 0.084 0.002 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 129 0.041 0.000 0.948 0.007 radau.py:49(solve_collocation_system)\n", | |
" 3004 0.040 0.000 0.040 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 1718 0.017 0.000 0.021 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 859 0.013 0.000 0.055 0.000 blas.py:226(_get_funcs)\n", | |
" 843 0.012 0.000 0.012 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 996 0.012 0.000 0.012 0.000 {method 'dot' of 'numpy.ndarray' objects}\n", | |
" 54 0.011 0.000 0.011 0.000 {built-in method numpy.core.multiarray.zeros}\n", | |
" 859 0.009 0.000 0.038 0.000 blas.py:177(find_best_blas_type)\n", | |
" 3716 0.009 0.000 0.016 0.000 shape_base.py:9(atleast_1d)\n", | |
" 358 0.006 0.000 0.013 0.000 linalg.py:1976(norm)\n", | |
" 1858 0.006 0.000 0.035 0.000 fromnumeric.py:1743(sum)\n", | |
" 981 0.005 0.000 0.005 0.000 {built-in method numpy.core.multiarray.arange}\n", | |
" 930 0.005 0.000 0.005 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
" 4073 0.005 0.000 0.006 0.000 numeric.py:484(asanyarray)\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10f6453c8>" | |
] | |
}, | |
"execution_count": 53, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## BDF method" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Small N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 54, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 3" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 55, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 20]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 56, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 58, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='BDF', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 59, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 710592 function calls in 1.298 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 91 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 5688 0.217 0.000 0.417 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 1960 0.116 0.000 1.279 0.001 bdf.py:286(_step_impl)\n", | |
" 2441 0.105 0.000 1.003 0.000 bdf.py:36(solve_bdf_system)\n", | |
" 25809 0.096 0.000 0.096 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 12448 0.092 0.000 0.121 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 6224 0.046 0.000 0.240 0.000 blas.py:226(_get_funcs)\n", | |
" 5686 0.046 0.000 0.322 0.000 decomp_lu.py:75(lu_solve)\n", | |
" 8410 0.045 0.000 0.091 0.000 linalg.py:1976(norm)\n", | |
" 22754 0.039 0.000 0.073 0.000 shape_base.py:9(atleast_1d)\n", | |
" 6224 0.031 0.000 0.183 0.000 blas.py:177(find_best_blas_type)\n", | |
" 12550 0.027 0.000 0.027 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 30402 0.024 0.000 0.030 0.000 numeric.py:484(asanyarray)\n", | |
" 1126 0.023 0.000 0.035 0.000 bdf.py:18(compute_R)\n", | |
" 13591 0.023 0.000 0.094 0.000 fromnumeric.py:1743(sum)\n", | |
" 227 0.022 0.000 0.025 0.000 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 79604 0.021 0.000 0.021 0.000 numerictypes.py:951(<listcomp>)\n", | |
" 6224 0.020 0.000 0.050 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 8410 0.019 0.000 0.110 0.000 common.py:38(norm)\n", | |
" 20551 0.018 0.000 0.023 0.000 numeric.py:414(asarray)\n", | |
" 5690 0.018 0.000 0.018 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10f6d1160>" | |
] | |
}, | |
"execution_count": 59, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Medium N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 65, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 20" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 66, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 20]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 67, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 68, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='BDF', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 69, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 816765 function calls in 1.893 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 91 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 6461 0.386 0.000 0.630 0.000 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 2200 0.146 0.000 1.864 0.001 bdf.py:286(_step_impl)\n", | |
" 2784 0.131 0.000 1.390 0.000 bdf.py:36(solve_bdf_system)\n", | |
" 29402 0.124 0.000 0.124 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 6459 0.119 0.000 0.437 0.000 decomp_lu.py:75(lu_solve)\n", | |
" 641 0.115 0.000 0.152 0.000 decomp_lu.py:17(lu_factor)\n", | |
" 14200 0.105 0.000 0.138 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 7100 0.054 0.000 0.277 0.000 blas.py:226(_get_funcs)\n", | |
" 9706 0.053 0.000 0.110 0.000 linalg.py:1976(norm)\n", | |
" 25846 0.043 0.000 0.082 0.000 shape_base.py:9(atleast_1d)\n", | |
" 245 0.037 0.000 0.041 0.000 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 14437 0.037 0.000 0.037 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 7100 0.036 0.000 0.209 0.000 blas.py:177(find_best_blas_type)\n", | |
" 7100 0.036 0.000 0.071 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 1442 0.030 0.000 0.045 0.000 bdf.py:18(compute_R)\n", | |
" 34507 0.028 0.000 0.034 0.000 numeric.py:484(asanyarray)\n", | |
" 15462 0.027 0.000 0.120 0.000 fromnumeric.py:1743(sum)\n", | |
" 90426 0.024 0.000 0.024 0.000 numerictypes.py:951(<listcomp>)\n", | |
" 9706 0.023 0.000 0.133 0.000 common.py:38(norm)\n", | |
" 6463 0.022 0.000 0.022 0.000 {built-in method numpy.core.multiarray.concatenate}\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10ee7fe10>" | |
] | |
}, | |
"execution_count": 69, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "markdown", | |
"metadata": {}, | |
"source": [ | |
"## Large N" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 76, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"N = 200" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 77, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"t_span=[0, 1]" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 78, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [ | |
"u0, m = setup(N)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 79, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" " | |
] | |
} | |
], | |
"source": [ | |
"stats = %prun -r -q solve_ivp(lambda t, y: fun(t, y, m), t_span, u0, method='BDF', jac=lambda t, y: jac(t, y, m))" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": 80, | |
"metadata": { | |
"collapsed": false | |
}, | |
"outputs": [ | |
{ | |
"name": "stdout", | |
"output_type": "stream", | |
"text": [ | |
" 32418 function calls in 1.829 seconds\n", | |
"\n", | |
" Ordered by: internal time\n", | |
" List reduced from 91 to 20 due to restriction <20>\n", | |
"\n", | |
" ncalls tottime percall cumtime percall filename:lineno(function)\n", | |
" 40 0.527 0.013 0.601 0.015 decomp_lu.py:17(lu_factor)\n", | |
" 258 0.471 0.002 0.500 0.002 <ipython-input-4-419434ad83ce>:1(fun)\n", | |
" 256 0.286 0.001 0.306 0.001 decomp_lu.py:75(lu_solve)\n", | |
" 74 0.239 0.003 1.804 0.024 bdf.py:286(_step_impl)\n", | |
" 25 0.098 0.004 0.127 0.005 <ipython-input-5-6e801fc71e23>:1(jac)\n", | |
" 296 0.070 0.000 0.075 0.000 function_base.py:970(asarray_chkfinite)\n", | |
" 76 0.026 0.000 0.026 0.000 {built-in method numpy.core.multiarray.zeros}\n", | |
" 1167 0.022 0.000 0.022 0.000 {method 'reduce' of 'numpy.ufunc' objects}\n", | |
" 108 0.015 0.000 0.834 0.008 bdf.py:36(solve_bdf_system)\n", | |
" 74 0.010 0.000 1.814 0.025 base.py:147(step)\n", | |
" 627 0.008 0.000 0.008 0.000 {built-in method numpy.core.multiarray.dot}\n", | |
" 592 0.006 0.000 0.007 0.000 numerictypes.py:942(_can_coerce_all)\n", | |
" 368 0.005 0.000 0.010 0.000 linalg.py:1976(norm)\n", | |
" 296 0.004 0.000 0.018 0.000 blas.py:226(_get_funcs)\n", | |
" 1034 0.003 0.000 0.005 0.000 shape_base.py:9(atleast_1d)\n", | |
" 296 0.003 0.000 0.012 0.000 blas.py:177(find_best_blas_type)\n", | |
" 600 0.003 0.000 0.021 0.000 fromnumeric.py:1743(sum)\n", | |
" 368 0.002 0.000 0.013 0.000 common.py:38(norm)\n", | |
" 399 0.002 0.000 0.002 0.000 {built-in method numpy.core.multiarray.arange}\n", | |
" 949 0.002 0.000 0.002 0.000 numeric.py:414(asarray)\n", | |
"\n", | |
"\n" | |
] | |
}, | |
{ | |
"data": { | |
"text/plain": [ | |
"<pstats.Stats at 0x10ee9bf60>" | |
] | |
}, | |
"execution_count": 80, | |
"metadata": {}, | |
"output_type": "execute_result" | |
} | |
], | |
"source": [ | |
"stats.print_stats(20)" | |
] | |
}, | |
{ | |
"cell_type": "code", | |
"execution_count": null, | |
"metadata": { | |
"collapsed": true | |
}, | |
"outputs": [], | |
"source": [] | |
} | |
], | |
"metadata": { | |
"anaconda-cloud": {}, | |
"kernelspec": { | |
"display_name": "Python [default]", | |
"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.5.2" | |
} | |
}, | |
"nbformat": 4, | |
"nbformat_minor": 1 | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment