{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "Contenido bajo licencia Creative Commons BY 4.0 y código bajo licencia MIT. © Juan Gómez y Nicolás Guarín-Zapata 2020. Este material es parte del curso Modelación Computacional en el programa de Ingeniería Civil de la Universidad EAFIT." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Interpolación de Lagrange en 1D" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Introducción" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "El problema de interpolación o de predecir funciones a partir de un número determinado de datos representativos de la función aparece con bastante frecuencia en la interpretación de datos experimentales. De otro lado, las técnicas de aproximación de funciones por medio de métodos de interpolación son la base para la formulación de los métodos más importantes de la mecánica computacional como lo son los elementos finitos y los elementos de frontera. En este Notebook se discutirán algunos aspectos básicos y fundamentales de la teoría de interpolación. El notebook se describe en términos del desarrollo completo de un problema ejemplo incluyendo su implementación en Python.\n", "\n", "**Al completar este notebook usted debería estar en la capacidad de:**\n", "\n", "* Reconocer el problema de interpolación de funciones como uno de aproximación de una función desconocida en términos de valores discretos de la misma función.\n", "\n", "* Identificar las propiedades fundamentales de los polinomios de interpolación de Lagrange.\n", "\n", "* Usar polinomios de interpolación de Lagrange para proponer aproximaciones a una función dado un conjunto de $N$ valores conocidos de la misma.\n", "\n", "* Reconocer la diferencia entre variables primarias y secundarias en un esquema de interpolación." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Interpolación de Lagrange\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "El problema de interpolación consiste en la determinación del valor de una función $f(x)$ en un punto arbitrario $x \\in [x_1, x_n]$, dados valores conocidos de la función al interior de un dominio de solución. De acuerdo con el teorema de interpolación de Lagrange la aproximación $\\hat f(x)$ a la función $f(x)$ se construye como:\n", "\n", "\\begin{equation}\n", "\\hat{f}(x) = \\sum_{I=1}^n L_I(x) f_I\n", "\\end{equation}\n", "\n", "donde $L_I$ es el $I-$ésimo polinomio de orden $n-1$ y $f_1, f_2, \\cdots, f_n$ son los $n$ valores conocidos de la función. El $I-$ésimo polinomio de Lagrange se calcula siguiendo la siguiente expresión:\n", "\n", "\\begin{equation}\n", "{L_I}(x) = \\prod_{J=1, J \\ne I}^{n} \\frac{(x - x_J)}{(x_I - x_J)}\\, .\n", "\\end{equation}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Primera derivada" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "En varios problemas de ingeniería es necesario usar técnicas de interpolación para encontrar valores de las derivadas de la variable primaria o principal. Por ejemplo, en problemas de mecánica de sólidos y teoría de elasticidad, en los cuales la variable primaria es el campo de desplazamientos, uno puede estar interesado en determinar las deformaciones unitarias las cuales son derivadas espaciales de los desplazamientos. Considerando que solo se dispone de valores discretos de los desplazamientos se tiene que las derivadas que se requieren se pueden encontrar usando estos valores discretos. Lo anterior equivale a derivar $\\hat{f}(x)$ directamente como sigue:\n", "\n", "\n", "\\begin{equation}\n", "\\frac{\\mathrm{d}\\hat{f}}{\\mathrm{d}x}=\\frac{\\mathrm{d}L_I(x)}{\\mathrm{d}x}f_I\n", "\\end{equation}\n", "\n", "Haciendo,\n", "\n", "$$B_I(x) = \\frac{\\mathrm{d}L_I(x)}{\\mathrm{d}x}\\, ,$$\n", "\n", "podemos escribir el esquema de interpolación como:\n", "\n", "\\begin{equation}\n", "\\frac{\\mathrm{d}\\hat{f}}{\\mathrm{d}x} = \\sum_{I=1}^n B_I(x) f_I\\, .\n", "\\end{equation}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Ejemplo\n", "\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Formule un esquema de interpolación para encontrar un valore de la función\n", "\n", "\\begin{equation}\n", "f(x) = x^3 + 4x^2 - 10\n", "\\end{equation}\n", "\n", "en un punto arbitrario $x$ en el intervalo $[ -1, 1]$ asumiendo que se conoce el valor exacto de la función en los puntos $x=-1.0$, $x=+1.0$ y $x=0.0.$\n", "\n", "En este ejemplo se conoce la función y aparentemente no tiene mucho sentido buscar una aproximación de la misma resolviendo un problema de interpolación. Sin embargo es conveniente seleccionar una función conocida para poder asimilar el problema numérico y sus limitaciones. En este contexto asumiremos que en una aplicación real conocemos los valores de la función en un conjunto de puntos $x=-1.0$, $x=+1.0$ and $x=0.0$ los cuales se denominan puntos nodales o simplemente _nodos_.\n", "\n", "El proceso de interpolación involucra 2 pasos fundamentales: \n", "\n", "1. Determinar los polinomios de interpolación $L_I$ usando la productoria.\n", "\n", "2. Usar la combinación lineal para construir el polinomio interpolante o la aproximación a la función $\\hat f(x)$.\n", "\n", "Veamos estos pasos entonces.\n", "\n", "1. Considerando que tenemos 3 puntos _nodales_ necesitamos generar 3 polinomios de interpolación de segundo orden, cada uno de ellos asociado a cada punto nodal. Rotulemos los _nodos_ como $x_0 = -1.0$, $x_1 = +1.0$ y $x_2 = 0.0$. De acuerdo con esta denominación $L_0(x)$, $L_1(x)$ y $L_2(x)$ serán los polinomios de interpolación de segundo orden asociados a los puntos nodales $x_0 = -1.0$, $x_1 = +1.0$ y $x_2 = 0.0$. Usando la fórmula de la productoria tenemos:\n", "\n", "\\begin{align}\n", "&L_0(x) = \\frac{(x - x_1)(x - x_2)}{(x_0 - x_1)(x_0 - x_2)} \\equiv \\frac{1}{2}(x - 1.0)x\\\\\n", "&L_1(x) = \\frac{(x-x^0)(x-x^2)}{(x^1-x^0)(x^1-x^2)}\\equiv\\frac12(x+1.0)x\\\\\n", "&L_2(x) = \\frac{(x-x^0)(x-x^1)}{(x^2-x^0)(x^2-x^1)}\\equiv-(x+1.0)(x-1.0)\\, .\n", "\\end{align}\n", "\n", "2. Para llegar a la aproximación final de la función realizamos la combinación lineal:\n", "\n", "\\begin{equation}\n", "\\hat{f}(x) = L_0(x)f_0 + L_1(x)f_1 + L_2(x)f_2\n", "\\end{equation}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "**Preguntas**\n", " \n", "- Verificar que los polinomios de interpolación $L_0(x)$, $L_1(x)$ y $L_2(x)$ satisfacen la propiedad\n", "\n", "$$L_I(x_J)= \\delta_{IJ} \\equiv \\begin{cases}\n", "1\\quad \\text{si } I=J\\\\\n", "0\\quad \\text{si } I\\neq J\n", "\\end{cases}\\, .$$\n", "\n", "- Si la condición $L_I(x_J)=\\delta_{IJ}$ no se satisface por una de las funciones de interpolación encontradas, ¿cuál sería el efecto resultante en la aproximación de la función?\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Algunas notas de interés\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "* En el método de los elementos finitos es común denominar a los polinomios de interpolación como _funciones de forma_.\n", "\n", "* El dominio de calculo es el intervalo de tamaño 2 comprendido entre $x=-1.0$ y $x=+1.0$. Como se discutirá mas adelante, por facilidades en la programación en la implementación de métodos de elementos finitos los dominios de tamaño diferente son transformados al dominio de tamaño 2." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solución en Python\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "En los siguientes bloques de código mostramos el paso a paso en la construcción del polinomio de interpolación final $f(x)$ usando Python. Para seguir el Notebook se recomienda que de manera simultánea se implementen los bloques de código en un script independiente o que añada comentarios a las instrucciones más relevantes." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Paso 1: Importación de módulos\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "En la escritura de scripts en Python es necesario importar **módulos** o bibliotecas (algunas personas usan la palabra _librería_ como mala traducción de _library_ del inglés) que contienen funciones de Python predefinidas. En este caso importaremos los **módulos**:\n", "\n", "* `numpy` el cual es una biblioteca de funciones para realizar operaciones con matrices similar a Matlab.\n", "* `matplotlib` el cual es una biblioteca de graficación.\n", "* `scipy` el cual es una biblioteca fundamental para computación científica.\n", "* `sympy` el cual es una biblioteca para realizar matemáticas simbólicas.\n", "\n", "Python importa los módulos usando la palabra reservada `import` seguida del nombre del módulo y un nombre corto para ser usado como prefijo en referencias posteriores a las funciones contenidas en ese módulo." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib notebook" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from scipy import interpolate\n", "import sympy as sym" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "sym.init_printing()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Paso 2: Definición de una función para determinar los polinomios de interpolación de Lagrange\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "En un programa de computador una **función** (o también llamada **subrutina**) es un bloque de código que realiza una tarea específica múltiples veces dentro de un programa o probablemente en diferentes programas. En Python estas funciones se definen por medio de la palabra clave `def` seguida del nombre dado a la función.\n", "\n", "Conjuntamente con el nombre, encerrado entre paréntesis, la definición de la función debe incluir también una lista de parámetros (o argumentos) a ser usados por la función cuando esta realice tareas especificas.\n", "\n", "En el ejemplo definiremos una función de Python para generar el polinomio de Lagrange usando la productoria definida previamente. Le daremos a esta función el nombre `lagrange_poly`. Sus parámetros de entrada serán la variable independiente $x$ a ser usada en el polinomio; el orden del polinomio definido por `order`; y el punto `i` asociado al polinomio. La función será usada posteriormente desde el programa principal." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "def lagrange_poly(x, order, i, xi=None): \n", " if xi == None:\n", " xi = sym.symbols('x:%d'%(order+1))\n", " index = list(range(order+1))\n", " index.pop(i)\n", " return sym.prod([(x - xi[j])/(xi[i] - xi[j]) for j in index])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Alternativamente a la definición de una función usando la palabra clave `def` Python también permite la definición de funciones en el sentido del Calculo, es decir definición de relaciones que permiten transformar escalares o vectores en otros escalares o vectores. Esto es posible a través de la opción `lambda`. En este ejemplo usaremos la opción `lambda` para crear la función:\n", "\n", "\\begin{equation}\n", "f(x)=x^3+4x^2-10.\n", "\\end{equation}\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "**Preguntas**\n", "\n", "Use una terminal o un script independiente para probar el uso de la opción `lambda` en la definición de una función. Intente con diferentes funciones.\n", "
" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [ "fun = lambda x: x**3 + 4.0*x**2 - 10.0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Una vez creada la función $f(x)$ usando la opción `lambda` pasamos a definir un conjunto de puntos de evaluación. El numero de puntos se define por medio de la variable `npts` y usamos la función `linspace` del módulo `numpy` para crear una distribución equidistante de puntos entre $x = -1.0$ y $x = +1.0.$ Simultáneamente, el arreglo nulo `yy` almacenará los valores de la función en los puntos almacenados en `xx`.\n", "\n", "Note que Python inicia desde la posición 0 el conteo de elementos en arreglos y otras estructuras de datos, de manera que empezamos a contar desde cero para mantener la consistencia con la implementación." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [], "source": [ "npts = 200\n", "x_pts = np.linspace(-1, 1, npts)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Con la función ahora disponible podemos calcular los valores (conocidos) listos para ser usados en el esquema de interpolación. Estos valores se almacenarán en el arreglo `fd()`. Para calcular cada valor de la función usamos el comando (ya disponible) `fx` correspondiente al nombre que usamos en la declaración de la función usando la opción `lambda`.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "**Preguntas**\n", "\n", "Intente usar nombres diferentes en la definición de la función usando la opción `lambda` y ensaye si su código aún funciona.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [], "source": [ "pts = np.array([-1, 1, 0])\n", "fd = fun(pts)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Ahora evalúe la función en los `npts` puntos y grafique la misma:" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support. ' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " if (mpl.ratio != 1) {\n", " fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n", " }\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " fig.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var backingStore = this.context.backingStorePixelRatio ||\n", "\tthis.context.webkitBackingStorePixelRatio ||\n", "\tthis.context.mozBackingStorePixelRatio ||\n", "\tthis.context.msBackingStorePixelRatio ||\n", "\tthis.context.oBackingStorePixelRatio ||\n", "\tthis.context.backingStorePixelRatio || 1;\n", "\n", " mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width * mpl.ratio);\n", " canvas.attr('height', height * mpl.ratio);\n", " canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
');\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.figure()\n", " \n", "for k in range(3):\n", " for i in range(npts):\n", " yy[i] = pol[k].subs([(x, xx[i])]) \n", " plt.plot(xx, yy)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "### Paso 4: Encontrando el polinomio de interpolación para aproximar la función $f(x)$\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Construyamos ahora el polinomio de aproximación completo $\\hat f(x)$ de acuerdo con la superposición lineal\n", "\n", "$$\\hat{f}(x) = L_0(x)f(x_0) + L_1(x) f(x_1) + L_2(x) f(x_2)\\, ,$$\n", "\n", "\n", "y utilizando la lista (ya disponible) `pol` que almacena los 3 polinomios generados. Solo para efectos de ilustrar primero lo mostramos con `display`. La versión evaluada de $\\hat{f}(x)$ será almacenada en el arreglo `yy[i]`." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAAZCAYAAAAyqZWXAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAIq0lEQVR4Ae2c7ZHWNhDHj5sr4CAVJHRAQgehAxIqIHQAwze+MdABpAICHQAV8NIBUAEMHVz+P+HVyD7bWj+2H78gzegk2bur1V/r1etzVy4uLk62Hh49evREbThX/E3xs+IDPfuutISCQEGgIFAQqBA42zoSlbN/phRHf6L0pZK3ir9TLqEgUBAoCBQEfiBwugMgbjfa8FjlG3L8zPZLKAgUBAoCBYEKgT04fJqSOnfbykmflQ4vCBQECgI/PQJ72NK53ujFG1X5feN5KRYECgIFgZ8agb3M8NNO5AC3HNqmiGTybH8p/pMh28xr2lK29A7rrmILh+G2Fa4rQ27pyBi4CYNDZdvkkyKz6xd6/lFpNozlz1Ug+eG2jtJ7Odq53qtutpLS+llxcKj8ylOn6D6I7oXic+W/K8KPvNdeGaJ1B8lEX/S75WbaAKHaw+E9A384zF9C5QrbB1XdtsWITt7vpdjCBB23BlvwNkO6zupja1s6qowD0H8Vf1Xe9sJTXbn98ljvgvOqlPug9Jai58May5/qUsurfmao50pTZ1ujmbuguvmo+aCjDsqD6UueKT536IAMBq4nojfyp8q7BgxjGJDiGP8aQL8KUuGRs9W7UnSx21rSL3y4SiO2ytOv9r28cQBZbMEB0tptIW2CQ9exPrKX/1QK4CRxSM+k2B1FDPVS0HscKkv/6HiUZ1CgjNPoDWP5+4RL9p96f11pcLSUFW0vv4916nc4IbYTSENQ3vCymZ696koZOJ8qwgcP7fLyitwfJJc+fa/UM1j7Bc9EKT1dtkr1osU23yi9P5M6ObFMnBh0YpAu9CN6Zb+XimnXtiA8bh/aP+LbjC14dRXd7D72TJVggGEWojwfR3RWldFZAk3brOSdnt+vGoWsrjCWv1Wu6sWxM3Ni5WG64/jjzKqVcZ6HLNXBIOIALlVV8Vmm6s/imcXBt9RLPZvZyhEuYOixVWvqY2W+KDKAHjswCWFLpnmpgG8IR8fV4dzWTrGFjl7bki0M0HWsj8zyn3bg2fYYA24L5si63htP13svv8lppixhcPrMmiwywze5TfrZyqqTGeVV0qQSazcrqNUE6ch2wYnSTczuDwFObcMGcJrWB4eIOZQHXFkRnx8q4Fh80rHYwrHA7q+ny07Nl3W9N6ld7yP/mVH2pYnRfuuhC0bT9n4of0XP8oZwU5GZKPJx7L8ofhVNmLUpvaryKkPVDlYf7MF79u9DO0RrbcdZMEPkULU2G6xkG10vRkHo5T+shtKBqUYxgfyavAULtLFr5jObWsKv65fettXoGmglx/q42ML43lrEFjxqV98bpLP62FOPMqK5VtHZSJGymYI44q7g5q8a/lApThKn/rqKLIEp46juKK42SE9mduGATkpiZGwteAMf9n/i55YO7bWDvjh66zk0YzFiK4cth0thIvmX5C704JPq/WOhumvVClf6kIkLtt32LdXoVSi20ERkXHk1ttDSDLePbOHlkYvf6/CtDgxwTPDwP1QFqYNkhstHYjNkZvvM2FYb9DGHvVel6M6H/aX62LM6i44bT/CEoDwzQQaNdEtoCozA1AbrUFfyZwr5ibhFs7TRY3fHUJI+fKU+xYazodhCFqKhBGuyhS7dx9pqL79rS0eadTkGlLaR5WtXCwbyc/gaHZ54mRWxNx6eKbVbLz3V+V5JFuBwBtALUiINHe6K72PyrDcrWq5psiznJlTXdddeGXqJ0+dcgpUD+Skwot8CpkqbYbR86Tk7tk2lO8rgxeDWG+bWV/Jx9tgxFwrGhC3aAqvUuEJNGo+NXBMmbSt2Jk1TT+xWYQtJ+9PsUXysy+EL+O+KKGfOPVWUTiMAZmsYwg9tQ8gNlV0zogZftljV1bXXmuVvEkhecCxKm1i8Fy0GT+wcsMTHoTNOvamTYYJ8PgQrqxjCpBhNIb+S0WyH6bu6dE59JTtcDVXqdvai3ZMttH6/aiPbs9h7OI9bi1FIH76vo9oudSoCwaw+9pQanIFthfMWWlOQ931hML8AwJERcrJ/UC3/lz3CT9K7DSePdgwIbcHkMXDUwgiMXMvbEfJrei5YwD6bA/DR1BF+OLX4GxEq1jOcXG7VUWxh+l5a1BYczRnsIxsys/xDHD4zjrbDL0bCjzJgRsW+4OKXHO4om4MLRq9y3EJRPtz576towXdgELefEj3s447tSN6lWQ5r22YW4BAxFs0UGOEEbbBOdTiZSH5N5oIFsM/Z5izqCUcmLJzJNGf2DAK5UGwhh9Dw94vZglNVl4/skZXldzt8GS2Hpt9wBlah8jjmvxXjLwp5pnihWLsBonKWXzQ4NpS22U3tR0F6HzpM6SIfsPTKBZau6B+DdAUv9I7/10XPWjESDf9Ogb3eGFRm/x/HHPYzVZ4KIwafS4PLhPJjGxbO0MajrxCFI32OLbBUp19D3yrl2T2lYdWhtNiCADlSWMQWvG2TLWR9JLK6bMbDH/55mghxMjgVnAlOHGeAQfIPu+x2zInyvHuoaAe0N5XngK82c1WZrQ1mpLVDlxx/9Z6fpL9TJFA3Do/bLmEAEU3UR89WF6QfGKbt5sPnY685HZW7MAqDg3gY1Mh/U2SwCIOcUvpgNEaSg57c7wfbGKaSHwVOnJF+Llu1akWP3Qw6aDfeMWlVLzP8tsA5TMRd+Z/SFtTuMBlSetAevvg2YQsYgFdX0c3qYwf9t8w2yy3PtouAjAtHw5bDYnvcc6KndjFgMmmJznXO+rYsewlbUJ2jHP4QvIst/EDrdAhohXZ3CHBdrvUGxU5ayt45s8AS8ggsYQusejtvreVVHkRRbEFwlRn+IJvZH7FmPvySOe4p76WFahdL47dKL51T7KWNU7ej2MLUiK5PXpnhr69Pjq3RXmc+HI6mZynHxnWL9RVb2GKvDdC5zPAHgLVXUs3s2Ovml7yrPhD34q92cNDP9dhdnk14cTiErtjCIahth+d/F7Qi5+fskz4AAAAASUVORK5CYII=\n", "text/latex": [ "$\\displaystyle 10.0 x^{2} - 3.5 x \\left(x - 1\\right) - 2.5 x \\left(x + 1\\right) - 10.0$" ], "text/plain": [ " 2 \n", "10.0⋅x - 3.5⋅x⋅(x - 1) - 2.5⋅x⋅(x + 1) - 10.0" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "display(pol[0]*fd[0] + pol[1]*fd[1] + pol[2]*fd[2])" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "application/javascript": [ "/* Put everything inside the global mpl namespace */\n", "window.mpl = {};\n", "\n", "\n", "mpl.get_websocket_type = function() {\n", " if (typeof(WebSocket) !== 'undefined') {\n", " return WebSocket;\n", " } else if (typeof(MozWebSocket) !== 'undefined') {\n", " return MozWebSocket;\n", " } else {\n", " alert('Your browser does not have WebSocket support. ' +\n", " 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n", " 'Firefox 4 and 5 are also supported but you ' +\n", " 'have to enable WebSockets in about:config.');\n", " };\n", "}\n", "\n", "mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n", " this.id = figure_id;\n", "\n", " this.ws = websocket;\n", "\n", " this.supports_binary = (this.ws.binaryType != undefined);\n", "\n", " if (!this.supports_binary) {\n", " var warnings = document.getElementById(\"mpl-warnings\");\n", " if (warnings) {\n", " warnings.style.display = 'block';\n", " warnings.textContent = (\n", " \"This browser does not support binary websocket messages. \" +\n", " \"Performance may be slow.\");\n", " }\n", " }\n", "\n", " this.imageObj = new Image();\n", "\n", " this.context = undefined;\n", " this.message = undefined;\n", " this.canvas = undefined;\n", " this.rubberband_canvas = undefined;\n", " this.rubberband_context = undefined;\n", " this.format_dropdown = undefined;\n", "\n", " this.image_mode = 'full';\n", "\n", " this.root = $('
');\n", " this._root_extra_style(this.root)\n", " this.root.attr('style', 'display: inline-block');\n", "\n", " $(parent_element).append(this.root);\n", "\n", " this._init_header(this);\n", " this._init_canvas(this);\n", " this._init_toolbar(this);\n", "\n", " var fig = this;\n", "\n", " this.waiting = false;\n", "\n", " this.ws.onopen = function () {\n", " fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n", " fig.send_message(\"send_image_mode\", {});\n", " if (mpl.ratio != 1) {\n", " fig.send_message(\"set_dpi_ratio\", {'dpi_ratio': mpl.ratio});\n", " }\n", " fig.send_message(\"refresh\", {});\n", " }\n", "\n", " this.imageObj.onload = function() {\n", " if (fig.image_mode == 'full') {\n", " // Full images could contain transparency (where diff images\n", " // almost always do), so we need to clear the canvas so that\n", " // there is no ghosting.\n", " fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n", " }\n", " fig.context.drawImage(fig.imageObj, 0, 0);\n", " };\n", "\n", " this.imageObj.onunload = function() {\n", " fig.ws.close();\n", " }\n", "\n", " this.ws.onmessage = this._make_on_message_function(this);\n", "\n", " this.ondownload = ondownload;\n", "}\n", "\n", "mpl.figure.prototype._init_header = function() {\n", " var titlebar = $(\n", " '
');\n", " var titletext = $(\n", " '
');\n", " titlebar.append(titletext)\n", " this.root.append(titlebar);\n", " this.header = titletext[0];\n", "}\n", "\n", "\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "\n", "mpl.figure.prototype._root_extra_style = function(canvas_div) {\n", "\n", "}\n", "\n", "mpl.figure.prototype._init_canvas = function() {\n", " var fig = this;\n", "\n", " var canvas_div = $('
');\n", "\n", " canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n", "\n", " function canvas_keyboard_event(event) {\n", " return fig.key_event(event, event['data']);\n", " }\n", "\n", " canvas_div.keydown('key_press', canvas_keyboard_event);\n", " canvas_div.keyup('key_release', canvas_keyboard_event);\n", " this.canvas_div = canvas_div\n", " this._canvas_extra_style(canvas_div)\n", " this.root.append(canvas_div);\n", "\n", " var canvas = $('');\n", " canvas.addClass('mpl-canvas');\n", " canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n", "\n", " this.canvas = canvas[0];\n", " this.context = canvas[0].getContext(\"2d\");\n", "\n", " var backingStore = this.context.backingStorePixelRatio ||\n", "\tthis.context.webkitBackingStorePixelRatio ||\n", "\tthis.context.mozBackingStorePixelRatio ||\n", "\tthis.context.msBackingStorePixelRatio ||\n", "\tthis.context.oBackingStorePixelRatio ||\n", "\tthis.context.backingStorePixelRatio || 1;\n", "\n", " mpl.ratio = (window.devicePixelRatio || 1) / backingStore;\n", "\n", " var rubberband = $('');\n", " rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n", "\n", " var pass_mouse_events = true;\n", "\n", " canvas_div.resizable({\n", " start: function(event, ui) {\n", " pass_mouse_events = false;\n", " },\n", " resize: function(event, ui) {\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " stop: function(event, ui) {\n", " pass_mouse_events = true;\n", " fig.request_resize(ui.size.width, ui.size.height);\n", " },\n", " });\n", "\n", " function mouse_event_fn(event) {\n", " if (pass_mouse_events)\n", " return fig.mouse_event(event, event['data']);\n", " }\n", "\n", " rubberband.mousedown('button_press', mouse_event_fn);\n", " rubberband.mouseup('button_release', mouse_event_fn);\n", " // Throttle sequential mouse events to 1 every 20ms.\n", " rubberband.mousemove('motion_notify', mouse_event_fn);\n", "\n", " rubberband.mouseenter('figure_enter', mouse_event_fn);\n", " rubberband.mouseleave('figure_leave', mouse_event_fn);\n", "\n", " canvas_div.on(\"wheel\", function (event) {\n", " event = event.originalEvent;\n", " event['data'] = 'scroll'\n", " if (event.deltaY < 0) {\n", " event.step = 1;\n", " } else {\n", " event.step = -1;\n", " }\n", " mouse_event_fn(event);\n", " });\n", "\n", " canvas_div.append(canvas);\n", " canvas_div.append(rubberband);\n", "\n", " this.rubberband = rubberband;\n", " this.rubberband_canvas = rubberband[0];\n", " this.rubberband_context = rubberband[0].getContext(\"2d\");\n", " this.rubberband_context.strokeStyle = \"#000000\";\n", "\n", " this._resize_canvas = function(width, height) {\n", " // Keep the size of the canvas, canvas container, and rubber band\n", " // canvas in synch.\n", " canvas_div.css('width', width)\n", " canvas_div.css('height', height)\n", "\n", " canvas.attr('width', width * mpl.ratio);\n", " canvas.attr('height', height * mpl.ratio);\n", " canvas.attr('style', 'width: ' + width + 'px; height: ' + height + 'px;');\n", "\n", " rubberband.attr('width', width);\n", " rubberband.attr('height', height);\n", " }\n", "\n", " // Set the figure to an initial 600x600px, this will subsequently be updated\n", " // upon first draw.\n", " this._resize_canvas(600, 600);\n", "\n", " // Disable right mouse context menu.\n", " $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n", " return false;\n", " });\n", "\n", " function set_focus () {\n", " canvas.focus();\n", " canvas_div.focus();\n", " }\n", "\n", " window.setTimeout(set_focus, 100);\n", "}\n", "\n", "mpl.figure.prototype._init_toolbar = function() {\n", " var fig = this;\n", "\n", " var nav_element = $('
');\n", " nav_element.attr('style', 'width: 100%');\n", " this.root.append(nav_element);\n", "\n", " // Define a callback function for later on.\n", " function toolbar_event(event) {\n", " return fig.toolbar_button_onclick(event['data']);\n", " }\n", " function toolbar_mouse_event(event) {\n", " return fig.toolbar_button_onmouseover(event['data']);\n", " }\n", "\n", " for(var toolbar_ind in mpl.toolbar_items) {\n", " var name = mpl.toolbar_items[toolbar_ind][0];\n", " var tooltip = mpl.toolbar_items[toolbar_ind][1];\n", " var image = mpl.toolbar_items[toolbar_ind][2];\n", " var method_name = mpl.toolbar_items[toolbar_ind][3];\n", "\n", " if (!name) {\n", " // put a spacer in here.\n", " continue;\n", " }\n", " var button = $('');\n", " button.click(method_name, toolbar_event);\n", " button.mouseover(tooltip, toolbar_mouse_event);\n", " nav_element.append(button);\n", " }\n", "\n", " // Add the status bar.\n", " var status_bar = $('');\n", " nav_element.append(status_bar);\n", " this.message = status_bar[0];\n", "\n", " // Add the close button to the window.\n", " var buttongrp = $('
');\n", " var button = $('');\n", " button.click(function (evt) { fig.handle_close(fig, {}); } );\n", " button.mouseover('Stop Interaction', toolbar_mouse_event);\n", " buttongrp.append(button);\n", " var titlebar = this.root.find($('.ui-dialog-titlebar'));\n", " titlebar.prepend(buttongrp);\n", "}\n", "\n", "mpl.figure.prototype._root_extra_style = function(el){\n", " var fig = this\n", " el.on(\"remove\", function(){\n", "\tfig.close_ws(fig, {});\n", " });\n", "}\n", "\n", "mpl.figure.prototype._canvas_extra_style = function(el){\n", " // this is important to make the div 'focusable\n", " el.attr('tabindex', 0)\n", " // reach out to IPython and tell the keyboard manager to turn it's self\n", " // off when our div gets focus\n", "\n", " // location in version 3\n", " if (IPython.notebook.keyboard_manager) {\n", " IPython.notebook.keyboard_manager.register_events(el);\n", " }\n", " else {\n", " // location in version 2\n", " IPython.keyboard_manager.register_events(el);\n", " }\n", "\n", "}\n", "\n", "mpl.figure.prototype._key_event_extra = function(event, name) {\n", " var manager = IPython.notebook.keyboard_manager;\n", " if (!manager)\n", " manager = IPython.keyboard_manager;\n", "\n", " // Check for shift+enter\n", " if (event.shiftKey && event.which == 13) {\n", " this.canvas_div.blur();\n", " event.shiftKey = false;\n", " // Send a \"J\" for go to next cell\n", " event.which = 74;\n", " event.keyCode = 74;\n", " manager.command_mode();\n", " manager.handle_keydown(event);\n", " }\n", "}\n", "\n", "mpl.figure.prototype.handle_save = function(fig, msg) {\n", " fig.ondownload(fig, null);\n", "}\n", "\n", "\n", "mpl.find_output_cell = function(html_output) {\n", " // Return the cell and output element which can be found *uniquely* in the notebook.\n", " // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n", " // IPython event is triggered only after the cells have been serialised, which for\n", " // our purposes (turning an active figure into a static one), is too late.\n", " var cells = IPython.notebook.get_cells();\n", " var ncells = cells.length;\n", " for (var i=0; i= 3 moved mimebundle to data attribute of output\n", " data = data.data;\n", " }\n", " if (data['text/html'] == html_output) {\n", " return [cell, data, j];\n", " }\n", " }\n", " }\n", " }\n", "}\n", "\n", "// Register the function which deals with the matplotlib target/channel.\n", "// The kernel may be null if the page has been refreshed.\n", "if (IPython.notebook.kernel != null) {\n", " IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n", "}\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "[]" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dpol = []\n", "dpol.append(sym.diff(pol[0],x)) \n", "dpol.append(sym.diff(pol[1],x))\n", "dpol.append(sym.diff(pol[2],x))\n", "display(dpol)\n", "display(dpol[0]*fd[0] + dpol[1]*fd[1] + dpol[2]*fd[2])\n", "\n", "plt.figure()\n", "yy= fdx(xx)\n", "plt.plot(xx, yy ,'r--')\n", "plt.plot([-1, 1, 0], fc, 'ko')\n", "\n", "for i in range(npts):\n", " yy[i] = fd[0]*dpol[0].subs([(x, xx[i])]) + fd[1]*dpol[1].subs([(x, xx[i])])\\\n", " + fd[2]*dpol[2].subs([(x, xx[i])]) \n", "plt.plot(xx, yy)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Glosario de términos\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Fórmula de productoria:** Expresión usada para calcular los polinomios de interpolación de Lagrange.\n", "\n", "**Función de forma:** Nombre dado a un polinomio de interpolación en el contexto del método de los elementos finitos.\n", "\n", "**Punto nodal:** (También nodo). Nombre dado a los puntos específicos donde se conoce el valor de una función y usado en la construcción del esquema de interpolación.\n", "\n", "**Subrutina:** Bloque de código independiente que realiza una tarea de cómputo determinada dentro de un programa principal.\n", "\n", "**Matriz de interpolación:** Arreglo unidimensional o bidimensional que almacena las funciones de forma en un esquema de interpolación determinado." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Actividad para la Clase" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "El propósito de esta actividad es familiarizar al estudiante con el uso del método de interpolación de Lagrange en un contexto particular de la ingeniería, como el método de los elementos finitos. El método de elementos finitos utiliza dominios de interpolación de tamaño constante (llamados **elementos**) permitiendo el uso de funciones de interpolación fijos y al mismo tiempo favoreciendo la automatización en programas de computador. En una aplicación típica de elementos finitos los valores nodales de una función (por ejemplo el desplazamiento) son determinados tras resolver un sistema de ecuaciones algebraicas. Estos valores nodales son posteriormente usados para encontrar el desplazamiento a lo largo de todo el dominio del problema, y también para realizar cálculos posteriores y obtener información adicional del problema físico.\n", "\n", "**Siga los pasos que se indican a continuación para implementar, en un notebook independiente, un esquema de interpolación típico del método de elementos finitos. El esquema de interpolación resultante será el correspondiente a un solo elemento**\n", "\n", "* Asumiendo que un dominio de interpolación constante se define como $x \\in [-1.0, +1.0]$ con puntos nodales correspondientes a $x_0 = -1.0$, $x_1 = +1.0$, $x_2 = 0.0$, $x_3 = -0.25$ y $x_4 = +0.25$ use la función `lagrange_poly()` para generar las funciones de interpolación asociadas a estos 5 puntos nodales. Presente los polinomios resultantes en una gráfica.\n", "\n", "* Utilice las funciones de interpolación encontradas en el paso anterior, para definir una matriz $[N(x)]$, que se denominará en adelante **Matriz de interpolación**, de manera que la operación de interpolar se pueda realizar usando operaciones matriciales como:\n", "\n", " $$\\hat{f}(x) = [N(x)] \\{F\\}\\, ,$$\n", "\n", " y en la cual $F$ es un vector que almacena los valores nodales de la función.\n", "\n", "* En su Notebook imprima la versión simbólica de la matriz de interpolación.\n", "\n", "* En su Notebook grafique una función $f(x)$ (representando el desplazamiento de una barra) y su versión aproximada $\\hat{f}(x)$ usando el esquema matricial desarrollado en el numeral anterior. El código también debe graficar la primera derivada de la función.\n", "\n", "* Para realizar la interpolación en el Notebook programe una función o subrutina de Python que reciba como parámetros de entrada las coordenadas $x$ donde se desea calcular una aproximación de la función y el vector de valores nodales de la función de desplazamientos y entregue después de su ejecución el valor interpolado de la función." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Formato del notebook" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "La siguiente celda cambia el formato del Notebook." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "text/plain": [ "" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from IPython.core.display import HTML\n", "def css_styling():\n", " styles = open('./nb_style.css', 'r').read()\n", " return HTML(styles)\n", "css_styling()" ] } ], "metadata": { "celltoolbar": "Raw Cell Format", "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.6.10" }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 1 }