Skip to content

Instantly share code, notes, and snippets.

@nmayorov
Last active November 3, 2016 22:49
Show Gist options
  • Save nmayorov/a16f3d9e90cfc7c86a96b1bc9e68b19f to your computer and use it in GitHub Desktop.
Save nmayorov/a16f3d9e90cfc7c86a96b1bc9e68b19f to your computer and use it in GitHub Desktop.
Profiling of ODE solvers
Display the source blob
Display the rendered blob
Raw
{
"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