{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Introduction\n", "\n", "Your goal is to assess the trade-offs of stability, accuracy, and run-time of a least two methods applied to two ODE systems. \n", "\n", "The first ODE system is two proteins degrading. The ODEs are:\n", "\n", "\\begin{eqnarray*}\n", "\\frac{dy_1}{dt}&=&-\\alpha y_1\\\\\n", "\\frac{dy_2}{dt}&=&-\\beta y_2\\\\\n", "\\end{eqnarray*}\n", "\n", "where $\\alpha=0.1$ and $\\beta=0.2$.\n", "\n", "The second ODE system is the Lotka-Volterra model of predator-prey populations and their interactions. You can check out the Wikipedia entry to learn more, but what you need to know here is that it is a 2-ODE system and that, with the right parameter values, it leads to oscillations. The ODEs are:\n", "\n", "\\begin{align*}\n", "\\frac{dx}{dt} &= \\alpha x -\\beta x y\\\\\n", "\\frac{dy}{dt} &= -\\gamma y + \\delta x y\\\\\n", "\\end{align*}\n", "\n", "where $\\alpha=0.25$, $\\beta=0.01$, $\\gamma=1$, and $\\delta=0.01$.\n", "\n", "\n", "## The Project\n", "\n", "Consider three methods -- the forward Euler method, the explicit trapezoidal method, and an adaptive step-size forward Euler method (that uses the trapezoidal method to estimate the local error).\n", "\n", "Compare the performance of each of your three methods on the degradation and predator-prey systems. In your analysis, you will want to consider both runtime (we want our answer as quickly as possible) and accuracy (how close is the method to the \"true\" solution?). For the first two methods, be sure to consider the effect of step-size on accuracy and speed. For the third, be sure to analyze the step-sizes chosen and the effect of error tolerance on the step size. The ideal graph plots accuracy versus runtime for different methods used to solve the same model (one line for each method, one step-size for each point on the line, where the step-size for the adaptive method is the average step-size chosen by the method).\n", "\n", "\n", "## Extensions\n", "\n", "This project will be more interesting if you tackle one or more of these extensions. Feel free to develop your own extension. And, as always, a particularly thorough analysis is smiled upon.\n", "\n", "\n", "\n", "## Notes\n", "\n", "### Timing\n", "To time your method, use the %timeit magic command. You use it either for a single line of code, e.g. \n", "\n", "\n", "%timeit (t,y) = forwardEuler( protein_degradation, np.arange(0,50,0.01), np.array([1,1]), params)\n", "\n", "\n", "or to time an entire cell, e.g.\n", "\n", "%%timeit\n", "(t,y) = forwardEuler( protein_degradation, np.arange(0,50,0.01), np.array([1,1]), params) \n", "\n", "\n", "If you use it in cell mode, then be sure to put just the solving code in the cell. Put any plotting code in a different cell. \n", "\n", "Regardless of which mode you use, you probably need to rerun the solving code before you do so. Stephanie has found that timeit doesn't actually do the assignment into (t,y) the way she was expecting it to. The safest thing to do is use a cell for just timing, and then a cell for re-running and plotting.\n", "\n", "### Accuracy\n", "To compute the accuracy of your method, compare its output to the output of scipy.integrate.solve_ivp (which will stand in for the \"true solution\"). Set its relative tolerance so that it is \"tight\" (i.e. you will tolerate very little error; rtol=1e-8). Error calculations should be performed by comparing the solutions to $y$ at the same time steps, but you will likely be comparing the true solution to outputs from methods using lots of different time steps. The easiest way to handle this is to call scipy.integrate.solve_ivp for each unique set of timesteps you have. For example, when you use the adaptive step size function, you do not know ahead of time which timesteps you will have, so run that first, and use the output to set the t_eval parameter for solve_ivp:\n", "\n", "(t,y,steps) = forwardEulerAdaptiveStepSize( protein_degradation, 0, 50, \n", " np.array([1,1]), params, 1e-4)\n", "sol = scipy.integrate.solve_ivp(fun=protein_degradation, t_span=(0,50), y0=(1,1), \n", " method='LSODA', t_eval=t, args=[params], rtol=1e-8, atol=1e-8 )\n", "\n", " \n", "To compute the error, use one of the methods we described in Lecture 12 (first, compute the error at each time point, then summarize it over time)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "Begin your code and report here. Organize your report in the manner that best allows you to answer all the questions. Integrate your extensions with your main project.\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Start code and write-up here." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "Acknowledgements\n", "\n", "*Replace this text with links to external resources and thanks to people you worked with.*" ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }