122 \frametitle{Outline} |
122 \frametitle{Outline} |
123 \tableofcontents |
123 \tableofcontents |
124 % \pausesections |
124 % \pausesections |
125 \end{frame} |
125 \end{frame} |
126 |
126 |
127 \section{Integration} |
127 \section{Interpolation} |
128 |
128 |
129 \subsection{Quadrature} |
129 \begin{frame}[fragile] |
|
130 \frametitle{Interpolation} |
|
131 \begin{itemize} |
|
132 \item Let us begin with interpolation |
|
133 \item Let's use the L and T arrays and interpolate this data to obtain data at new points |
|
134 \end{itemize} |
|
135 \begin{lstlisting} |
|
136 In []: L = [] |
|
137 In []: T = [] |
|
138 In []: for line in open('pendulum.txt'): |
|
139 l, t = line.split() |
|
140 L.append(float(l)) |
|
141 T.append(float(t)) |
|
142 In []: L = array(L) |
|
143 In []: T = array(T) |
|
144 In []: Tsq = T*T |
|
145 \end{lstlisting} |
|
146 \end{frame} |
|
147 |
|
148 %% \begin{frame}[fragile] |
|
149 %% \frametitle{Interpolation \ldots} |
|
150 %% \begin{small} |
|
151 %% \typ{In []: from scipy.interpolate import interp1d} |
|
152 %% \end{small} |
|
153 %% \begin{itemize} |
|
154 %% \item The \typ{interp1d} function returns a function |
|
155 %% \begin{lstlisting} |
|
156 %% In []: f = interp1d(L, T) |
|
157 %% \end{lstlisting} |
|
158 %% \item Functions can be assigned to variables |
|
159 %% \item This function interpolates between known data values to obtain unknown |
|
160 %% \end{itemize} |
|
161 %% \end{frame} |
|
162 |
|
163 %% \begin{frame}[fragile] |
|
164 %% \frametitle{Interpolation \ldots} |
|
165 %% \begin{lstlisting} |
|
166 %% In []: Ln = arange(0.1,0.99,0.005) |
|
167 %% # Interpolating! |
|
168 %% # The new values in range of old data |
|
169 %% In []: plot(L, T, 'o', Ln, f(Ln), '-') |
|
170 %% In []: f = interp1d(L, T, kind='cubic') |
|
171 %% # When kind not specified, it's linear |
|
172 %% # Others are ... |
|
173 %% # 'nearest', 'zero', |
|
174 %% # 'slinear', 'quadratic' |
|
175 %% \end{lstlisting} |
|
176 %% \end{frame} |
|
177 |
|
178 \begin{frame}[fragile] |
|
179 \frametitle{Spline Interpolation} |
|
180 \begin{small} |
|
181 \begin{lstlisting} |
|
182 In []: from scipy.interpolate import splrep |
|
183 In []: from scipy.interpolate import splev |
|
184 \end{lstlisting} |
|
185 \end{small} |
|
186 \begin{itemize} |
|
187 \item Involves two steps |
|
188 \begin{enumerate} |
|
189 \item Find out the spline curve, coefficients |
|
190 \item Evaluate the spline at new points |
|
191 \end{enumerate} |
|
192 \end{itemize} |
|
193 \end{frame} |
|
194 |
|
195 \begin{frame}[fragile] |
|
196 \frametitle{\typ{splrep}} |
|
197 To find the B-spline representation |
|
198 \begin{lstlisting} |
|
199 In []: tck = splrep(L, T) |
|
200 \end{lstlisting} |
|
201 Returns a tuple containing |
|
202 \begin{enumerate} |
|
203 \item the vector of knots, |
|
204 \item the B-spline coefficients |
|
205 \item the degree of the spline (default=3) |
|
206 \end{enumerate} |
|
207 \end{frame} |
|
208 |
|
209 \begin{frame}[fragile] |
|
210 \frametitle{\typ{splev}} |
|
211 To Evaluate a B-spline and it's derivatives |
|
212 \begin{lstlisting} |
|
213 In []: Lnew = arange(0.1,1,0.005) |
|
214 In []: Tnew = splev(Lnew, tck) |
|
215 |
|
216 #To obtain derivatives of the spline |
|
217 #use der=1, 2,.. for 1st, 2nd,.. order |
|
218 In []: Tnew = splev(Lnew, tck, der=1) |
|
219 \end{lstlisting} |
|
220 \end{frame} |
|
221 |
|
222 \begin{frame}[fragile] |
|
223 \frametitle{} |
|
224 \end{frame} |
|
225 |
|
226 \section{Differentiation} |
|
227 |
|
228 \begin{frame}[fragile] |
|
229 \frametitle{Numerical Differentiation} |
|
230 \begin{itemize} |
|
231 \item Given function $f(x)$ or data points $y=f(x)$ |
|
232 \item We wish to calculate $f^{'}(x)$ at points $x$ |
|
233 \item Taylor series - finite difference approximations |
|
234 \end{itemize} |
|
235 \begin{center} |
|
236 \begin{tabular}{l l} |
|
237 $f(x+h)=f(x)+h.f^{'}(x)$ &Forward \\ |
|
238 $f(x-h)=f(x)-h.f^{'}(x)$ &Backward |
|
239 \end{tabular} |
|
240 \end{center} |
|
241 \end{frame} |
|
242 |
|
243 \begin{frame}[fragile] |
|
244 \frametitle{Forward Difference} |
|
245 \begin{lstlisting} |
|
246 In []: x = linspace(0, 2*pi, 100) |
|
247 In []: y = sin(x) |
|
248 In []: deltax = x[1] - x[0] |
|
249 \end{lstlisting} |
|
250 Obtain the finite forward difference of y |
|
251 \end{frame} |
|
252 |
|
253 \begin{frame}[fragile] |
|
254 \frametitle{Forward Difference \ldots} |
|
255 \begin{lstlisting} |
|
256 In []: fD = (y[1:] - y[:-1]) / deltax |
|
257 In []: plot(x, y, x[:-1], fD) |
|
258 \end{lstlisting} |
|
259 \begin{center} |
|
260 \includegraphics[height=2in, interpolate=true]{data/fwdDiff} |
|
261 \end{center} |
|
262 \end{frame} |
|
263 |
|
264 \begin{frame}[fragile] |
|
265 \frametitle{} |
|
266 \end{frame} |
|
267 |
|
268 \section{Quadrature} |
130 |
269 |
131 \begin{frame}[fragile] |
270 \begin{frame}[fragile] |
132 \frametitle{Quadrature} |
271 \frametitle{Quadrature} |
133 \begin{itemize} |
272 \begin{itemize} |
134 \item We wish to find area under a curve |
273 \item We wish to find area under a curve |
135 \item Area under $(sin(x) + x^2)$ in $(0,1)$ |
274 \item Area under $(sin(x) + x^2)$ in $(0,1)$ |
136 \item scipy has functions to do that |
275 \item scipy has functions to do that |
137 \end{itemize} |
276 \end{itemize} |
138 \small{\typ{In []: from scipy.integrate import quad}} |
277 \begin{small} |
|
278 \typ{In []: from scipy.integrate import quad} |
|
279 \end{small} |
139 \begin{itemize} |
280 \begin{itemize} |
140 \item Inputs - function to integrate, limits |
281 \item Inputs - function to integrate, limits |
141 \end{itemize} |
282 \end{itemize} |
142 \begin{lstlisting} |
283 \begin{lstlisting} |
143 In []: x = 0 |
284 In []: x = 0 |
144 In []: quad(sin(x)+x**2, 0, 1) |
285 In []: quad(sin(x)+x**2, 0, 1) |
145 \end{lstlisting} |
286 \end{lstlisting} |
|
287 \begin{small} |
146 \alert{\typ{error:}} |
288 \alert{\typ{error:}} |
147 \typ{First argument must be a callable function.} |
289 \typ{First argument must be a callable function.} |
|
290 \end{small} |
148 \end{frame} |
291 \end{frame} |
149 |
292 |
150 \begin{frame}[fragile] |
293 \begin{frame}[fragile] |
151 \frametitle{Functions - Definition} |
294 \frametitle{Functions - Definition} |
|
295 We have been using them all along. Now let's see how to define them. |
152 \begin{lstlisting} |
296 \begin{lstlisting} |
153 In []: def f(x): |
297 In []: def f(x): |
154 return sin(x)+x**2 |
298 return sin(x)+x**2 |
155 In []: quad(f, 0, 1) |
299 In []: quad(f, 0, 1) |
156 \end{lstlisting} |
300 \end{lstlisting} |
231 \item Use \typ{dblquad} for Double integrals |
375 \item Use \typ{dblquad} for Double integrals |
232 \item Use \typ{tplquad} for Triple integrals |
376 \item Use \typ{tplquad} for Triple integrals |
233 \end{itemize} |
377 \end{itemize} |
234 \end{frame} |
378 \end{frame} |
235 |
379 |
236 \subsection{ODEs} |
|
237 |
|
238 \begin{frame}[fragile] |
|
239 \frametitle{ODE Integration} |
|
240 We shall use the simple ODE of a simple pendulum. |
|
241 \begin{equation*} |
|
242 \ddot{\theta} = -\frac{g}{L}sin(\theta) |
|
243 \end{equation*} |
|
244 \begin{itemize} |
|
245 \item This equation can be written as a system of two first order ODEs |
|
246 \end{itemize} |
|
247 \begin{align} |
|
248 \dot{\theta} &= \omega \\ |
|
249 \dot{\omega} &= -\frac{g}{L}sin(\theta) \\ |
|
250 \text{At}\ t &= 0 : \nonumber \\ |
|
251 \theta = \theta_0\quad & \&\quad \omega = 0 \nonumber |
|
252 \end{align} |
|
253 \end{frame} |
|
254 |
|
255 \begin{frame}[fragile] |
|
256 \frametitle{Solving ODEs using SciPy} |
|
257 \begin{itemize} |
|
258 \item We use the \typ{odeint} function from scipy to do the integration |
|
259 \item Define a function as below |
|
260 \end{itemize} |
|
261 \begin{lstlisting} |
|
262 In []: def pend_int(unknown, t, p): |
|
263 .... theta, omega = unknown |
|
264 .... g, L = p |
|
265 .... f=[omega, -(g/L)*sin(theta)] |
|
266 .... return f |
|
267 .... |
|
268 \end{lstlisting} |
|
269 \end{frame} |
|
270 |
|
271 \begin{frame}[fragile] |
|
272 \frametitle{Solving ODEs using SciPy \ldots} |
|
273 \begin{itemize} |
|
274 \item \typ{t} is the time variable \\ |
|
275 \item \typ{p} has the constants \\ |
|
276 \item \typ{initial} has the initial values |
|
277 \end{itemize} |
|
278 \begin{lstlisting} |
|
279 In []: t = linspace(0, 10, 101) |
|
280 In []: p=(-9.81, 0.2) |
|
281 In []: initial = [10*2*pi/360, 0] |
|
282 \end{lstlisting} |
|
283 \end{frame} |
|
284 |
|
285 \begin{frame}[fragile] |
|
286 \frametitle{Solving ODEs using SciPy \ldots} |
|
287 |
|
288 \small{\typ{In []: from scipy.integrate import odeint}} |
|
289 \begin{lstlisting} |
|
290 In []: pend_sol = odeint(pend_int, |
|
291 initial,t, |
|
292 args=(p,)) |
|
293 \end{lstlisting} |
|
294 \end{frame} |
|
295 |
|
296 \begin{frame} |
380 \begin{frame} |
297 \frametitle{Things we have learned} |
381 \frametitle{Things we have learned} |
298 \begin{itemize} |
382 \begin{itemize} |
|
383 \item Interpolation |
|
384 \item Differentiation |
299 \item Functions |
385 \item Functions |
300 \begin{itemize} |
386 \begin{itemize} |
301 \item Definition |
387 \item Definition |
302 \item Calling |
388 \item Calling |
303 \item Default Arguments |
389 \item Default Arguments |
304 \item Keyword Arguments |
390 \item Keyword Arguments |
305 \end{itemize} |
391 \end{itemize} |
306 \item Quadrature |
392 \item Quadrature |
307 \end{itemize} |
393 \end{itemize} |
308 \end{frame} |
394 \end{frame} |
309 |
|
310 \end{document} |
395 \end{document} |
311 |
396 |