{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Intégration numérique en Python\n", "\n", "##### Germain Salvato-Vallverdu [germain.vallverdu@univ-pau.fr](germain.vallverdu@univ-pau.fr)\n", "\n", "L'objectif de ce TP est de mettre en pratique le langage python pour intégrer numériquement une fonction.\n", "Après avoir étudier la fonction des méthodes basiques comme la méthode des trapèzes ou celle de Simpson\n", "seront implémentées. Puis les fonctions disponnibles dans numpy et scipy seront utilisées pour réaliser\n", "les mêmes opérations. La librairie matplotlib sera utilisée tout au long du notebook pour faire des\n", "représentations graphiques.\n", "\n", "Note : Ce notebook est compatible python2 et python3. Il est recommandé d'utiliser python3.\n", "\n", "### Sommaire\n", "\n", "1. [Modules python](#Modules-python)\n", "2. [Étude d'une fonction](#Étude-d'une-fonction)\n", "3. [Intégration analytique](#Intégration-analytique)\n", "4. [Intégration numérique](#Intégration-numérique)\n", " 1. [Rappels mathématiques](#Rappels-mathématiques)\n", " 2. [À la main](#À-la-main)\n", " 3. [Avec numpy](#Avec-numpy)\n", "\n", "### Modules python\n", "\n", "Commençons par charger les modules nécessaires : matplotlib pour les représentations graphiques et numpy pour les approches numériques. Ces modules peuvent être chargés selon\n", "\n", "```python\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "```\n", "Avec ipython ou jupyter cela peut se faire via une commande magique. Ces commandes commence par un %. L'option `inline` permet d'intégrer les graphiques dans le notebook." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "%pylab --no-import-all inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note sur la documentation :\n", "\n", "```\n", " --no-import-all Prevent IPython from performing `import *` into the interactive namespace.\n", "```" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Étude d'une fonction\n", "\n", "Dans cet exemple, nous étudierons la fonction suivante :\n", "\\begin{equation}\n", " f(x) = (x^2 - 3x - 6) \\exp\\left(-x\\right)\n", "\\end{equation}\n", "Cette fonction est intégrable analytiquement. Il est alors nécessaires de faire deux intégrations par partie successives.\n", "\n", "Commençons par définir une fonction qui renvoie l'image de la fonction $f$ pour une valeur de $x$. La fonction exponnentielle n'est pas disponnible dans python par défaut, il faut la prendre dans le module numpy :" ] }, { "cell_type": "code", "execution_count": 112, "metadata": {}, "outputs": [], "source": [ "def fonction(x):\n", " \"\"\" Retourne la valeur de la fonction en x \"\"\"\n", " return (x**2 - 3*x - 6) * np.exp(-x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Utilisons maintenant matplotlib pour représenter la fonction entre dans l'intervalle $[0;5]$." ] }, { "cell_type": "code", "execution_count": 113, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 113, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAENCAYAAAAc1VI3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHqpJREFUeJzt3Xl0VFW2x/EvIraoSETRVkTppagtqMEZcAjIpCD61EZR\nsOM8oTI4YHejtMpgQzcO2E+cxQGeoiAyC1gIAipDbGcFiZHZVhBklKTeH7sSSqwklcqtOvfe+n3W\nqmXdqpvKdq9a7Jyz7z0HREREREREREREREREREREREREREREsl4Nh7+7DfA/wFogCty/y/t7AkOB\n5UBjYDDwdSYDFBER9/bC/vGvFTseA7Te5Zy+wB2x502BdzMTmoiIJLKbo9/bHPgW+CV2/B7QcZdz\nzgPmxZ5/ApwA7JOR6ERE5DdcFYwDgY1xxxtir1X1HBERyRBXBWMNUCfuuG7stXhrgX3jjveNvSYi\nIg7s7uj3zgcOB/YAtgMtgMeB/YAd2MhiIjZ1NQc4DigAft71gw455JDoypUrMxO1iEg4LAWOrOoP\nub5K6hLge6xoPAA8BPwY+2/pVVKrsP+xAcCSBJ8TjUajmYg39Pr370///v1dhxEayqe3lE/v1KhR\nA1L499/VCANgeuwR7+6451uBHpkLRwoLC12HECrKp7eUT/dc9TBERCRgVDCkTH5+vusQQkX59Jby\n6Z7LHoZX1MMQEamCVHsYGmFImUgk4jqEUFE+vaV8uqeCISIiSdGUlIhIltGUlIiIpJUKhpTRHLG3\nlE9vKZ/uqWCIiEhS1MMQEcky6mGIiEhaqWBIGc0Re0v59Jby6Z4KhoiIJEU9DBGRLKMehoiIpJUK\nhpTRHLG3lE9vKZ/uqWCIiEhS1MMQEckyQdyiVUREMmTVKli4EBYsSP0zNCUlZTRH7C3l01vKZ/LW\nrYNp02DAALjwQmjQAJo2hcceg+3bU/9cjTBERAJs2zYoKID334cPPrDH6tVw4olwyilw+eUwbBg0\nagQ1YpNQgwal9rvUwxARCYhoFIqKYN48mD/fHh9/DI0bw+mnw6mnwmmnwTHHQM2a5X9Oqj0MFQwR\nEZ/atg0WL4a5c3c+SkqgeXMrEKefDiefDHvvXbXPVcGQaotEIuTl5bkOIzSUT29lQz7XrbOiMGeO\nPRYvttFDixY7H/FTS6nSVVIiIgGzYgXMng3vvmv/LSy0KaWWLeHee+35vvu6jnInjTBERDKksBBm\nzdr5+OknOOMMOOssOPNMyM2FWrXSH4empEREfObbbyES2fnYsgXOPtseZ50Fxx4Luzm4uUEFQ6ot\nG+aIM0n59FYQ8rlmDcycufOxcSO0amWPvDw4+ujq9x+8oB6GiEiGbdxoU0vTp8OMGbB8uY0eWreG\n22+HJk38USC84uJ/pR4wCPgGaAz8BVib4LxCYFns+XKgezmfpxGGiGTEjh22tMa0afYoKLDGdJs2\ncM45drPc7gH4MzxIU1JPANOBMUAnoAtwZYLz7gP+nsTnqWCISNoUFcHUqfaYORMaNoR27aBtW2tY\n77WX6wirLkgFowhoDqzARhtfA/snOG86MBWoA0wG5pXzeSoYHgnCHHGQKJ/eylQ+t22zy1wnT4Yp\nU+D77604tG9vheLgg9MeQtr5rYcxBTgowev3AgcCG2PHG4D9sEUQS3Y5ty+wAKgNLMJGI0vTEayI\nZLeiIpg0CSZOtJ5E06Zw7rnwwgtw0klurmTyIz+PMOKNwkYZIxO8pxGGiFRJcbGtwzRhgj1Wr4YO\nHaBjRxtN7F/Zv0gB57cRRkUmAi2A14CWwITY6zWAhlhBaQ3UwqakAI4ElpT3gfn5+TRq1AiAnJwc\ncnNzy4aupUsi61jHOs7u4w0b4F//ijB3LixenEeDBnDccRFuvBFuvDGPmjXt/I8/9ke8Xh6XPi8s\nLKQ6XIww9gMeAr4FjgDuBr4HcrERxPFAU6A/sBA4BBuNDC7n8zTC8EhEc+6eUj69lUo+v/sOxo+3\nx7x5tuTG+edDp05w2GHpiTMIgjTCWAdcn+D1AqxYAHwCXJKxiEQkFKJR+OwzGDsWxo2DZcvgvPPg\nuutgzBioU8d1hMEWhltKNMIQyWLRqN0b8frr8MYbsHWr7TJ34YW2PlMm1mYKmiBdVus1FQyRLFNS\nYsuAjxljRaJ2bbj4YrjoIruqKUx3V6dDqgVDF4tJmfgGmVSf8umtGTMizJoFPXrAoYfCLbdAvXp2\nv8QXX8DAgbaZkIpF+gTgJnYRyVYlJfDee/Dqq/Dyy3D44fCnP9nKr0cd5Tq67BOGWqwpKZEQKe1J\njB4N//d/sN9+0KWLPY4+2nV04RCkq6RERH7j88/hlVdg1CibVura1dZvatLEdWRSSj0MKaM5d28p\nn5VbsQKGDoVmzWzF182bbWTx1Vdw//2/LhbKp3saYYhIRm3caJfAvvQSLFxoVzb985+2j0TNmq6j\nk4qohyEiaVdcbJsMvfCCLfB39tnQvbvdcV27tuvoso/uwxAR3/n8c3j+eRtNHHww/PnPcNllUL++\n68iym+7DkGrTHLG3sjWf69fDE0/YTnTnnGNXPU2bZlc+3Xpr6sUiW/PpJ+phiEi1lZTYPhLPPGPL\nhbdtC/fdZxsOBWHLUkmOpqREJGUrV8Jzz8Gzz9pWpddcA926wQEHuI5MKqL7MEQkI4qLbevSJ5+E\n2bPtzutRo+CUU7QsR9iphyFlNEfsrbDlc/lyuzfiD3+w/3bubFubjhgBp56a/mIRtnwGkUYYIlKu\nkhK7HPZ//9d6FJdeapsR5ea6jkxcCMMAUj0MEY/98IP1Jp54wjYduukmW6pDGxCFg3oYIlJtCxbA\n44/bjnWdO8OLL8Lpp6s3IUY9DCmjOWJvBSWf27fbon/Nm9smRMccA19/DSNH2mt+KRZByWeYaYQh\nkqXWrLEppxEj4I9/hLvvhvPP13pOUj6f/O1QLephiFTB4sXw8MPWvO7SBW67TUuIZxstDSIi5Sou\nhnHjbNG/zp3h2GNh6VIbXahYSLJUMKSM5oi95Yd8btoEw4fbTnWDBtnVTt98Y9NP9eq5jq5q/JDP\nbKcehkgIrV5thWLECDjzTFtWvEUL/zSwJZjC8PVRD0Mk5ssvbQe711+3ZcR794Yjj3QdlfiN7sMQ\nyWLz58PgwTBvHtx8s21xqgUAxWvqYUgZzRF7K935jEZh8mRrZHftakuKL1tmy4qHsVjo++meRhgi\nAVNcbFNOAwfaWk99+9rlsdp3QtJNPQyRgPjlF9vqdNAg2H9/+OtfoWNHNbKl6tTDEAmprVttIcCH\nHoLGje3Kp7w8FQrJPBc9jBrADcAa4NgKzusGDAUeAq7PQFxZT3PE3qpuPrdsgUcftaucJk2C0aPh\n7behVavsLBb6frrnYoRxAjAf2FzBOYcCfYBmseMPgJnAkvSGJuLe5s02ihgyxDYmGj8eTjzRdVQi\nbnsYy4COwGcJ3rsGaA5cGzt+BCsWjyU4Vz0MCYUtW2wxwH/8w1aJvfdebVQk6eG3HsYU4KAEr98L\nvJXEz9cHNsYdbwAO9CAuEd/ZutX2xx482PaemDxZhUL8KV0Fo0M1f34tEH9/al3gq/JOzs/Pp1Gj\nRgDk5OSQm5tLXl4esHPeU8eVH8fPEfshnqAfV5bP7duhb98IL74IzZvnMXEi/PRThPXrAdzH77dj\nfT9TPy59XlhYSHW4npLqBHwaF0tDoAhoAEzg1z2MrsDSBJ+jKSmPRCKRsi+aVF95+dyxwy6P/fvf\n4aij4IEHrFchFdP30zupTkm5KBg5QA+gF/AiMAp4H8gFRgLHx867AjgZKAa+BJ4q5/NUMCQQSkrg\njTegXz+oXx8GDLCFAUUyLUgFw2sqGOJr0ShMnw733GPPBw6Edu2y89JY8QdtoCTVFj/fKdUXiUT4\n4ANo0wZ69LAlPBYsgPbtVSxSoe+ne7rTWyQNvv7aFgFcutT+e9VVWutJgi8Mf+doSkp8Y80aa2a/\n+ir06QO33w577eU6KpFf05SUiEObNtnVTsceC3vsAV98YT0LFQsJExUMKaM54qorLoZnn7U9sz/9\nFD78EB5+2PajUD69pXy6p1lVkRTNmGHTTvvsY/tTnHaa64hE0ks9DJEq+uoruOMOG1H84x9w0UW6\n6kmCRT0MkTRbtw5694aWLeGss+Czz+Dii1UsJHuoYEgZzREnVlxsq8gec4wtPf7ppzbC+N3vKv45\n5dNbyqd76mGIVCASsUtj69WzzYuOP77SHxEJrTAMptXDEM8VFVlD+8MPYehQTT1JuKiHIeKBrVvh\nwQehWTNo2hQ+/xwuuUTFQgRUMCROts8RT5wITZrAokWwcKEt6VG7duqfl+359Jry6Z56GJL1li2D\nnj1tNPHvf9vigCLyW2EYaKuHISnZtg2GDIFhw6xf0adP5Vc+iYSB3/b0FvG1GTPg5pvtUtmFCyG2\nw6+IVEA9DCmTDXPEq1fD5ZfDNdfY6OLNN9NXLLIhn5mkfLqngiFZoaQERoyA446Dhg3t5rvOnV1H\nJRIs6mFI6H3yCdxwg22PWlo0RLKZ7sMQ2cXWrfC3v0GrVtC9O8yZo2IhUh0qGFImTHPE77xjy3h8\n+SV89BHceCPsluFve5jy6QfKp3u6SkpCZf16uPNOmDIFhg+HCy5wHZFIeKiHIaExbhz06GHN7EGD\noG5d1xGJ+JPuw5CstXYt3HorFBTAK6/YXhUi4j31MKRM0OaIo1ErEMcfb/dSFBT4q1gELZ9+p3y6\npxGGBNKqVXDTTbBkCbz1FpxyiuuIRMIv2Tms44AmwJ5AEfAusCNdQVWRehhZpHRU0bs3XHcd9Oun\n9Z9EqirVHkZlP9AAeAErFCuAX4B6wOHA3cCEqv7CNFDByBKrV9vlsUuXwvPPw0knuY5IJJjSceNe\nDaAP0B04A7gU6AacB5wANAOOreovFP/y8xzxa69Bbq7tV7FgQTCKhZ/zGUTKp3sV9TCiWMFI9Of7\n3sADpNY0rwFcD9wPtAI+K+e8QmBZ7PlyrHBJlvnhB7jlFmtov/kmnHaa64hEsleyQ5J+WIEo9Sxw\ndYq/MxcrQuOAjpRfMO4D/p7E52lKKqQmT4Zrr4UuXWDgwOrtficiO6X7Poz2wAZgJPA8Nh2VqoIk\nzzsTuBOoA0wG5lXjd0qA/Pwz3HGHFYyXXrK1oETEvWSnlM4HfsIa35OxnkZFpgCLEzzOr0JsfYEh\nwCBsRHNEFX5WUuCHOeL586FZM1s48D//CXax8EM+w0T5dC/ZEcY72GW0bYBLgKOBXhWc36GacQEs\niP13CzYqaQksTXRifn4+jWK74OTk5JCbm0teXh6w80umY38fn3FGHg8+CI88EqFnT7jvPn/Fp2Md\nB/m49HlhYSHVkewc1iTgImBr7GeeBa6q1m+2hnYn4NO4WBpi93m0BmoBU2PvfQjcDsxN8DnqYQTc\nkiXQrRvk5MCzz8Ihh7iOSCTc0r0fxk1YsYCdV0+lKgf4G7AvcB1Qet3LCey8r2Nt7L17gMeA10lc\nLCTAolG7n6J5c9s2ddIkFQsRP6uowuwGDIg9fk7w/i3YHd8fpyGuqtAIwyORSKRsKJtu69bZLnif\nf253bodxY6NM5jMbKJ/eScdVUiXYvRAzsWb3CqyPUQ8bDQzDfbGQAJozB664wvaqGDkS9tzTdUQi\nkozKKswIbDqqE/AHYA/gO+ArYFF6Q0uaRhgBsWMHDBgATzwBTz8NHTu6jkgkO6XrPoxtwGHYHdnD\n4n7BLfinYEgAfPedjSr22AMWLYKDD3YdkYhUVWVN7/nAk8CF2CKEz8cel6Q1KnEi/hI8L40fDyef\nDOedB9OmZU+xSFc+s5Xy6V5lI4xXYo/OwPi41zWZIJXatg3uusvWgBo7Flq0cB2RiFSH9vSWtFi6\nFC69FBo2hGeegXr1XEckIqXSfR+GSNLGjLF7K668Et54Q8VCJCxUMKRMdeeIt22DW2+1aaiJE+G2\n26BGGMawKdKcu7eUT/e0p7d4orDQliFv0MCugsrJcR2RiHgtDH//qYfh2IQJcM01cPfd0KtXdo8q\nRIIg3fthiPzGjh1w773w4ou6CkokG6iHIWWqMke8Zg20awcffmhTUCoWv6U5d28pn+6pYEiVzZ0L\nJ50ELVvClClQv77riEQkE8Iw26weRoZEo/D443D//bZvRadOriMSkVSohyFptXmzLUf+n//AvHlw\nhDbMFck6mpKSMuXNES9bZtNP0aiKRVVozt1byqd7KhhSobfftru28/Ptaqi99nIdkYi4oh6GJBSN\nwpAhMGwYjBoF2uhMJDzUwxDPbN5sN+J9/TW8/z4cdpjriETEDzQlJWUikQjffmv9it13h9mzVSyq\nQ3Pu3lI+3VPBkDIffQSnnw7du9te27Vru45IRPxEPQwBYMQIW+bjpZegbVvX0YhIOqmHISn55Rfo\n2RNmzoQ5c6BxY9cRiYhfaUoqi/3wA7Rvb0uTz58PK1ZEXIcUKppz95by6Z4KRpb64gvrV5x4Iowf\nD3Xruo5IRPxOPYwsNG0adOsGgwfD1Ve7jkZEMk09DEnKv/9tiweOGQNnneU6GhEJEk1JZYkdO2yP\n7cces+XJExULzRF7S/n0lvLpnkYYWWDDBrjsMisa8+Zpv20RSY2LHsYwYBPwM3AC0BNYk+C8bkAu\nUAwsBZ4s5/PUw6hAUZHtW9GihY0uatVyHZGIuBakHsbPQL/Y87uAvwK37XLOoUAfoFns+ANgJrAk\nEwGGxcKF0Lkz9OkDvXpBjTBc4iAizrjoYfSLe14T2JjgnPbAwrjjecC56QwqbN58Ezp0gOHDoXfv\n5IqF5oi9pXx6S/l0L10jjCnAQQle7wdMiD3PAdoCFyU4rz6/LiQbgAO9DDDMHnsMBg2CSZPglFNc\nRyMiYZGugtGhkvfrAsOBq4D1Cd5fCxy5y/lflfdh+fn5NGrUCICcnBxyc3PJi23gUPpXSTYcl5TA\nZZdFmD8f5s7No1Gjqv18Xl6er/5/gn6sfCqffjkufV5YWEh1uJjVPgB4GOtfrAQuBl6PxdIQKAIa\nYCOR+B5GV6z5vSs1vYEtW2yV2f/+F8aOhf32cx2RiPhVqk1vFz2MqUBT4GXgHeDK2OsnsHO6agUw\nFLuiaijwFImLhWBrQrVpA3vsAVOnpl4s4v8akepTPr2lfLrn4iqpk8p5vQA4Pu745dhDKrBsmTW3\nL7zQ+ha76VZMEUmTMFxombVTUqWXzd5zD/To4ToaEQmKIN2HIR6YNg2uuMI2Proo0XVmIiIe0wRG\nAL38sjW4x471tlhojthbyqe3lE/3NMIImH/+Ex591HbIa9LEdTQikk3UwwiIkhK46y6YMsUehx7q\nOiIRCSr1MELsl1/gmmtg6VJ4912oV891RCKSjdTD8LlNm+yS2R9/hLffTm+x0Byxt5RPbymf7qlg\n+Ni6ddCuHRxwgDW499rLdUQiks3Uw/CpVausWLRtC0OH6oY8EfFOkJYGkUosXQpnnAFdu9pVUSoW\nIuIH+qfIZz75BM4+G+68E/7yl8xueqQ5Ym8pn95SPt3TVVI+8v77ttTHI4/YHtwiIn6iHoZPzJxp\nReK556BjR9fRiEiY6T6MAHvrLbvP4rXXbDpKRMSP1MNwbPRouO46mDjRfbHQHLG3lE9vKZ/uqWA4\n9Mwz0KcPTJ+uvbdFxP/Uw3Dk0Uftktnp06FxY9fRiEg2UQ8jQAYPhqeftnWhDj/cdTQiIsnRlFQG\nRaPQvz+88ALMmuW/YqE5Ym8pn95SPt3TCCNDolHbSnXSJCsWBx7oOiIRkapRDyMDolHo1Qtmz7at\nVfff33VEIpLN1MPwqZIS6NEDFi2CGTMgJ8d1RCIiqVEPI41KSuCGG+Cjj2xk4fdioTlibymf3lI+\n3dMII02Ki+Haa23l2SlToE4d1xGJiFSPehhpUFwMV18NRUUwYQLsvbfriEREdlIPwyeKiyE/H1au\ntOU+tEueiISFehgeKi6GK6+03fLeeit4xUJzxN5SPr2lfLqnguGR0mKxdi2MHx+8YiEiUhn1MDxQ\nOg21apWKhYj4X5B6GMOATcDPwAlAT2BNgvMKgWWx58uB7pkIrqpKG9wrVwZzGkpEJFkupqR+Bv4G\nDAYWA38t57zngFaxhy+LRUmJXTr73XfhKBaaI/aW8ukt5dM9FyOMfnHPawIbyznvTOBOoA4wGZiX\n5riqpPSmvG++sfWhgl4sREQqk64exhTgoASv9wMmxJ7nAG8AFwHrE5x7MrAAqA0sAjoBSxOcl/Ee\nRjQKN98MH39sN+Xts09Gf72ISLX4rYfRoZL36wLDgatIXCzAigXAFqAAaEnigpFR0SjcfjsUFNhy\nHyoWIpItXExJHQA8DNwFrAQuBl7Hql1DoAhoDdQCpsZ+5khgSXkfmJ+fT6NGjQDIyckhNzeXvLw8\nYOe8pxfH0ShcdlmEggL44IM86tTx9vNdH8fPEfshnqAfK5/Kp1+OS58XFhZSHS4uq12I9S7WxY43\nABcAucBI4HigKdA/du4hwAqsSZ5Ixqak+vWz5vbMmVCvXkZ+ZUZFIpGyL5pUn/LpLeXTO6lOSek+\njCQNGACvvAKRCNSvn/ZfJyKSNn7rYYTKsGEwcqTtlKdiISLZSkuDVOKJJ+DRR2H6dPj9711Hk17x\n851Sfcqnt5RP9zTCqMDIkTYVNWsWNGzoOhoREbfUwyjHmDFw663wzjtwzDGef7yIiDPqYXho8mS4\n5RaYOlXFQkSklHoYu5g1y5YpHzcOcnNdR5NZmiP2lvLpLeXTPRWMOAsWwJ/+BKNHQ/PmrqMREfEX\n9TBiPvsMWreGJ5+Ezp09iEpExKdS7WFohAEsWwbt28PQoSoWIiLlyfqCsWoVtG0LfftCt26uo3FL\nc8TeUj69pXy6l9UFY906G1nk59tVUSIiUr6s7WFs2gTt2sGpp8K//gU1wpAJEZEkaPHBKti+HS64\nAA48EJ57DnbL6nGWiGQbNb2TVFJiU1C1asHTT6tYxNMcsbeUT28pn+5l1Z3e0Sj07AnLl9td3LVq\nuY5IRCQ4smpKasAAePVVu5s7JyfNUYmI+JTWkqrEU0/BM8/Ae++pWIiIpCIrZvDHjYP77rNpqIMP\ndh2Nf2mO2FvKp7eUT/dCP8KYPRuuvx6mTIHGjV1HIyISXKHuYXzyCZxzDrz8MrRpk+GoRER8SpfV\n7qKoCM49Fx5+WMVCRMQLoSwYP/4IHTpA797QtavraIJDc8TeUj69pXy6F7qCsWWLrTh73nnQq5fr\naEREwiNUPYziYujSBX73O3jpJd3FLSKSSNbfh1F6F/e6dbYnt4qFiIi3QvPP6pAhdgf32LE2wpCq\n0xyxt5RPbymf7oVihDFqFAwfDnPnQt26rqMREQmnUPQw6tePMmMGHHec61BERPwvq+/DGD1axUJE\nJN1cFIzbgaeAO4FxwOnlnNcNGAo8BFxf0Qe2bu1leNlLc8TeUj69pXy656Jg7AH0AIYAzwP3Jzjn\nUKAPcAdwN3AtcGSG4staBQUFrkMIFeXTW8qney4KxhBgW+x5Y+DTBOe0BxbGHc8Dzk1zXFlv/fr1\nrkMIFeXTW8qne+m6SmoKcFCC1/sBE2Lv/QXIBS5KcF59YGPc8QbgQI9jFBGRKkhXwehQyftrsF5G\nK2AScNou76/l11NQdYGvPItOEiosLHQdQqgon95SPrPTHXHP/wB8H3teAzgs9rwBsDjuvA+AI8r5\nvCVAVA899NBDj6QfS0iBi/swHgW2A/8FTgBGAeOx6amRwPGx864ATgaKgS+xK6tERERERERERCTw\ngrI0yO+BB7HpqlMTvL8ndpPfcuxS3cHA1xmLLlgqy2U+cAOwNXb8DPBSRiILpiOAB4BF2P1DP8SO\n4+n7mbxk8pmPvqPJqIFdlTofu//tCOBqduYNQvrdvBjoBHxYzvt92dlMbwq8m4mgAqqyXP4ZODxz\n4QTeycD5ccefAifuco6+n8lLJp/6jianBnb7QqlxwOW7nBPa72Ye5f8j9y7QMu74J2CfdAcUYHlU\nXDCGYHfa9wP2y1BMYfE5cNQur+n7mbpE+dR3tOp2Bxbw2+Jbpe9mKJY3x27qS3Sj389uwgm0Wdgw\n9gfs7vrXgDZOIwqO/8FuWt31niF9P1NTXj71Ha2adkAv4C1sqi9eaL+beVQ8wjgj7lh/wVUsj/Jz\nGW9P4BeC0+tyqRXwcDnv6ftZdRXlM56+o8l7Abhpl9eq9N0M8vLm+wF1Ys8nAs1jz48DCghBhcyg\n+FwOBGrGnjcGlmE3+kj5OmJ/xfUEDsFWYNb3M3WV5VPf0eT8ETgv7rgQu1k65e9mzfLe8JmzgO7Y\nzX21sb+O+2FNmvewhQq7AM2wL9sdwI9OIvW/8nLZBMtlE+xKiqZYg/wuYIWTSIPhJGx6pAZ29c6V\nwDfYdIq+n1WXTD71HU3Ovthq342Bc4BjgHuxrSX03RQRERERERERERERERERERERERERERERERER\nERGRcLkH2AycCfTGll5o7DQiEY8EZWkQkaCYgy1rcSq2iNv1wFqnEYmIiG/VwPZxuNl1ICIi4m8X\nY5v8fIqtDioiIvIbVwPvY5vQvIAtF32m04hERERERERERERERERERERERERERERERERERERExH/+\nH7GBW8pqu6HJAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "x = np.linspace(1, 3, 100)\n", "plt.plot(x, fonction(x))\n", "plt.grid()\n", "plt.xlabel(\"x\")\n", "plt.ylabel(\"f(x)\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Intégration analytique\n", "\n", "L'intégration analytique donne :\n", "\n", "\\begin{equation}\n", " \\int_1^3 (x^2 - 3x - 6) \\exp\\left(-x\\right) dx = -2.525369\n", "\\end{equation}" ] }, { "cell_type": "code", "execution_count": 83, "metadata": {}, "outputs": [], "source": [ "analytique = -2.525369" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Intégration numérique\n", "\n", "### Rappels mathématiques\n", "\n", "\n", "\n", "### Par la méthode des trapèzes\n", "\n", "Dans la méthode des trapèzes, l'intervalle d'intégration est découpé en petits intervalles. On dit que la fonction est discrétisée. Sur chaque petit intervalle, la fonction est approximée par une droite. L'intégrale est évaluée en calculant l'aire sous la courbe comme la somme des petites aires $S_i$, sur chaque petit intervalle, entre la droite qui approxime la courbe et l'axe des abscisses.\n", "\n", "![trapeze](./trapeze.png)\n", "\n", "Chaque point de la courbe est repéré par son indice $i$.\n", "\n", "L'aire d'un trapèze sur le petit intervalle compris entre $x_{i+1}$ et $x_i$ est donnée par :\n", "\\begin{equation}\n", " S_i = (x_{i+1} - x_i)\\frac{f(x_{i+1}) + f(x_i)}{2}\n", "\\end{equation}\n", "L'intégrale est alors la somme des petites aires $S_i$.\n", "\\begin{equation}\n", " I = \\sum_{i=1}^{N-1} S_i = \\sum_{i=1}^{N-1} (x_{i+1} - x_i)\\frac{f(x_{i+1}) + f(x_i)}{2}\n", "\\end{equation}\n", "Ce qui peut se simplifier de la façon suivante :\n", "\\begin{equation}\n", " I = h\\left(\\frac{f(a)+f(b)}{2} + \\sum_{i=2}^{N-1}f(x_i)\\right)\n", "\\end{equation}\n", "Avec $h$ la largeur des petits intervalles : le pas d'intégration." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Mise en œuvre de la méthode des trapèze." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "On définit les bornes d'intégration :" ] }, { "cell_type": "code", "execution_count": 84, "metadata": {}, "outputs": [], "source": [ "a, b = 1, 3" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "On choisit le nombre de valeurs de x :" ] }, { "cell_type": "code", "execution_count": 141, "metadata": {}, "outputs": [], "source": [ "n = 201" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calcul de $h$, le pas. Avec $n$ valeurs de $x$, on a $n-1$ intervalles." ] }, { "cell_type": "code", "execution_count": 142, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.01\n" ] } ], "source": [ "h = (b - a) / (n - 1)\n", "print(h)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "On calcule les valeurs des abscisses $x_i$ :" ] }, { "cell_type": "code", "execution_count": 154, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "201 1.0 1.01 3.0\n" ] } ], "source": [ "x = np.linspace(a, b, n)\n", "print(len(x), x[0], x[1], x[-1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Première possibilité :" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Voici une première façon de faire. Il faut bien penser à initialiser la variable contenant le résultat de l'intégration." ] }, { "cell_type": "code", "execution_count": 150, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-2.52538674545\n" ] } ], "source": [ "integrale = 0\n", "for i in range(n-1):\n", " integrale += (x[i+1] - x[i])*(fonction(x[i+1]) + fonction(x[i])) / 2\n", "print(integrale)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Seconde possibilité :" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Seconde façon de faire avec la formule simplifiée." ] }, { "cell_type": "code", "execution_count": 155, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-2.52538674545\n" ] } ], "source": [ "integrale = (fonction(a) + fonction(b)) / 2\n", "for xi in x[1:n-1]:\n", " integrale += fonction(xi)\n", "integrale *= h\n", "print(integrale)" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "ename": "NameError", "evalue": "name 'fonction' is not defined", "output_type": "error", "traceback": [ "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[1;31mNameError\u001b[0m Traceback (most recent call last)", "\u001b[1;32m\u001b[0m in \u001b[0;36m\u001b[1;34m()\u001b[0m\n\u001b[1;32m----> 1\u001b[1;33m \u001b[0mintegrale\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;33m(\u001b[0m\u001b[0mfonction\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0ma\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;33m+\u001b[0m \u001b[0mfonction\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mb\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;33m/\u001b[0m \u001b[1;36m2\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m 2\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0mxi\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m:\u001b[0m\u001b[0mn\u001b[0m\u001b[1;33m-\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 3\u001b[0m \u001b[0mintegrale\u001b[0m \u001b[1;33m+=\u001b[0m \u001b[0mfonction\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mxi\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 4\u001b[0m \u001b[0mintegrale\u001b[0m \u001b[1;33m*=\u001b[0m \u001b[0mh\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 5\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mintegrale\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n", "\u001b[1;31mNameError\u001b[0m: name 'fonction' is not defined" ] } ], "source": [ "integrale = (fonction(a) + fonction(b)) / 2\n", "for xi in x[1:n-1]:\n", " integrale += fonction(xi)\n", "integrale *= h\n", "print(integrale)\n" ] } ], "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": 1 }