diff --git a/Calcul_Diff/Doc/.DS_Store b/Calcul_Diff/Doc/.DS_Store new file mode 100644 index 0000000000000000000000000000000000000000..5008ddfcf53c02e82d7eee2e57c38e5672ef89f6 Binary files /dev/null and b/Calcul_Diff/Doc/.DS_Store differ diff --git a/Calcul_Diff/Doc/Cours-fin.pdf b/Calcul_Diff/Doc/Cours-fin.pdf new file mode 100644 index 0000000000000000000000000000000000000000..0d89d06b6595df5aabce434b82f047a3879d825e Binary files /dev/null and b/Calcul_Diff/Doc/Cours-fin.pdf differ diff --git a/Calcul_Diff/Doc/CoursLM256.pdf b/Calcul_Diff/Doc/CoursLM256.pdf new file mode 100644 index 0000000000000000000000000000000000000000..fc07d746f868c5450b68250b2047e8b3ff1f6150 Binary files /dev/null and b/Calcul_Diff/Doc/CoursLM256.pdf differ diff --git a/Calcul_Diff/Doc/int-curv.pdf b/Calcul_Diff/Doc/int-curv.pdf new file mode 100644 index 0000000000000000000000000000000000000000..2c2c12fa04e95577ce08c35d9729d78c93c39ca6 Binary files /dev/null and b/Calcul_Diff/Doc/int-curv.pdf differ diff --git a/Calcul_Diff/slides_calcul_diff.pdf b/Calcul_Diff/slides_calcul_diff.pdf new file mode 100644 index 0000000000000000000000000000000000000000..075018910a893b41dbc23a922781c3305db92e91 Binary files /dev/null and b/Calcul_Diff/slides_calcul_diff.pdf differ diff --git a/ode/TP/ode-raide-etudiant.ipynb b/ode/TP/ode-raide-etudiant.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..92404b6f96c3307e521473db4eb32c8bf97ae219 --- /dev/null +++ b/ode/TP/ode-raide-etudiant.ipynb @@ -0,0 +1,203 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Problème raide (stiff)" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "Plots.PyPlotBackend()" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# import Pkg; Pkg.add(\"DifferentialEquations\")\n", + "# Pkg.add(\"Plots\")\n", + "# import Pkg; Pkg.add(\"ODEInterfaceDiffEq\")\n", + "using DifferentialEquations\n", + "using Plots\n", + "using ODEInterfaceDiffEq # pour accéder à radau5\n", + "pyplot()\n", + "# gr()\n", + "# gr(dpi=300,size=(600,400),thickness_scaling=1)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "\n", + "## Exemple curtiss et Hirschfelder\n", + "Cet exemple provient du tome 2 du livre d'Hairer, page 2. Le problème de Cauchy est\n", + "$$(IVP)_2\\left\\{\\begin{array}{l}\n", + "\\dot{y}(t)=-50(y(t)-\\cos(t))\\\\\n", + "y(0)=0\n", + "\\end{array}\\right.\n", + "$$\n", + "avec $[t_0\\;\\; t_f]=[0\\;\\; 1.5]$.\n", + "\n", + "On demande de résoudre numériquement ce problème avec les algorithmes : \n", + "- d'Euler explicite, alg = Euler(), et avec un pas de 1.974/50 et un pas de 1.875/50 (option dt). On prendra aussi comme option dense = false afin de n'avoir que les points aux instant $t_i$.\n", + "- d'Euler implicite, algo = ImplicitEuler(), avec un pas de 0.5. On prendra aussi comme options dense = false et adaptive = false afin de n'avoir que les points aux instant $t_i$.\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "\u001b[36mODEProblem\u001b[0m with uType \u001b[36mArray{Float64,1}\u001b[0m and tType \u001b[36mFloat64\u001b[0m. In-place: \u001b[36mfalse\u001b[0m\n", + "timespan: (0.0, 1.5)\n", + "u0: [0.0]" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function curtiss(x,p,t)\n", + "# deuxieme membre de l'equation differentiel de l'equation de Curtiss et Hirschfelder\n", + "# ref: Hairer page 2 tome 2\n", + "# Input\n", + "# x : state\n", + "# real(2)\n", + "# p : parameter vector\n", + "# t : time\n", + "# real\n", + "# Output\n", + "# xpoint : vector of velocity\n", + "# same as x\n", + " xpoint = similar(x)\n", + " \n", + "# A compléter\n", + " return xpoint\n", + "end\n", + "\n", + "x0 = [0.];\n", + "t0 = 0.;\n", + "tf = 1.5;\n", + "p = [1.]\n", + "tspan = (t0,tf)\n", + "prob = ODEProblem(curtiss,x0,tspan,p)\n", + "# A compléter\n", + "# ..." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Exemple de Roberston\n", + "On considère la réaction chimique \n", + "$$\\begin{array}{cccl}\n", + "A & \\stackrel{0.04}{\\longrightarrow} & B & \\textrm{(lente)},\\\\\n", + "B+B & \\stackrel{3.10^7}{\\longrightarrow} & C+B & \\textrm{(très rapide)},\\\\\n", + "B+C & \\stackrel{10^4}{\\longrightarrow} & A+C & \\textrm{(rapide)}.\\\\\n", + "\\end{array}\n", + "$$\n", + "Le système différentiel associé à cette réaction chimique est donnée par\n", + "$$(IVP)\\left\\{\\begin{array}{lrr}\n", + "\\dot{x}_1(t)= &-0.04x_1(t)+10^4x_2(t)x_3(t) &\\\\\n", + "\\dot{x}_2(t)= & 0.04x_1(t)-10^4x_2(t)x_3(t) &-3.10^7x_2^2(t)\\\\\n", + "\\dot{x}_3(t)= & & 3.10^7x_2^2(t)\\\\\n", + "x_1(0)=1,\n", + "x_2(0)=0,\n", + "x_3(0)=0.\n", + "\\end{array}\\right.\n", + "$$\n", + "avec $[t_0\\;\\; t_f]=[0\\;\\; 0.3]$.\n", + "\n", + "On demande de résoudre numériquement ce problème avec les algorithmes : \n", + "- algo = DP5() (Dormand-Prince's 5/4 Runge-Kutta method, algorithme utilisé dans ode45 de matlab et dopri5 du Professeur Hairer) et les options reltol = 1e-2, abstol = 1.e-6, dense = false\n", + "- algo = DP5() (Dormand-Prince's 5/4 Runge-Kutta method, algorithme utilisé dans ode45 de matlab et dopri5 du Professeur Hairer) et les options reltol = 1e-3, abstol = 1.e-6, dense = false, dense = false\n", + "- algo = radau5() : pragramme Radau5 du Professeur Hairer pour les problème raides.\n", + "\n", + "Et d'afficher sur le même graphique la deuxième composante\n", + "\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "\u001b[36mODEProblem\u001b[0m with uType \u001b[36mArray{Float64,1}\u001b[0m and tType \u001b[36mFloat64\u001b[0m. In-place: \u001b[36mfalse\u001b[0m\n", + "timespan: (0.0, 0.3)\n", + "u0: [1.0, 0.0, 0.0]" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function chimie(x,p,t)\n", + "# deuxieme membre de l'equation differentiel de l'equation de Robertson (1966)\n", + "# ref: Hairer page 3 tome 2\n", + "# Input\n", + "# x : state\n", + "# real(2)\n", + "# p : parameter vector\n", + "# t : time\n", + "# real\n", + "# Output\n", + "# xpoint : vector of velocity\n", + "# same as x\n", + " xpoint = similar(x)\n", + "# A compléter\n", + " \n", + " return xpoint\n", + "end\n", + "\n", + "x0 = [1., 0., 0.];\n", + "t0 = 0.;\n", + "tf = 0.3;\n", + "p = []\n", + "tspan = (t0,tf)\n", + "xpoint = similar(x0)\n", + "prob = ODEProblem(chimie,x0,tspan,p)\n", + "# A compléter\n", + "# ..." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.4.1", + "language": "julia", + "name": "julia-1.4" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.4.1" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/ode/TP/ode-rk-etudiant.ipynb b/ode/TP/ode-rk-etudiant.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..1d8dd6cac2b057cd14343a147a1b46a275f88094 --- /dev/null +++ b/ode/TP/ode-rk-etudiant.ipynb @@ -0,0 +1,705 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# TP Algorithme de Runge Kutta\n", + "\n", + "**Document à rendre au plus tard le 25 septembre à 18h00**\n", + "\n", + "Nom : \n", + "\n", + "Prénom :\n", + "\n", + "## Introduction\n", + "\n", + "Il est rappelé que les programmes doivent respecter les **bonnes pratiques de la programmation**. En particulier on vérifiera que les interfaces soient bien définies (paramètres en entrée, en sortie avec leurs types, les dimensions, ...). Dans le cas contraire on mettra des **points négatifs** pour un maximum de 4 points.\n", + "\n", + "On rappelle qu'un schéma de Runge Kutta à $s$ étage est défini par\n", + "\\begin{equation}\\label{eq:irk}\n", + "(S)\\left\\{\\begin{array}{l}\n", + "k_i=f(t_0+c_ih,x_0+h\\sum_{j=1}^{s}a_{ij}k_j)\\quad \\textrm{pour}\\;\\;i=1,\\ldots,s\\\\\n", + "x_1=x_0+h\\sum_{i=1}^{s}b_ik_i\n", + "\\end{array}\\right.\n", + "\\end{equation}\n", + "où les coefficients $c_i,a_{ij}$ et $b_i$ sont des constantes qui définissent précisément le schéma. On supposera toujours dans la suite que $c_i=\\sum_{j=1}^{s}a_{ij}$ pour $i=1,\\ldots,s$. \n", + "\n", + "On représente en pratique ce schéma par le tableau de Butcher\n", + "$$\n", + "\\begin{array}{c|ccc}\n", + "c_1 & a_{11} & \\ldots & a_{1s}\\\\\n", + "\\vdots & \\vdots & & \\vdots\\\\ \n", + "c_s & a_{s1} & \\ldots & a_{ss} \\\\ \\hline\n", + " & b_1 & \\ldots & b_s\\\\\n", + "\\end{array}$$\n", + "\n", + "Si la matrice des coefficient $A$ est triangulaire inférieure stricte alors le schéma est explicite, sinon il est implicite.\n", + "\n", + "On donne ci après les schémas les plus courant\n", + "\n", + "$$\\begin{array}[t]{ccc}\n", + "\\begin{array}{c}\n", + "\\\\ \\\\\n", + "\\begin{array}{c|c}\n", + "0 & \\\\\\hline\n", + " & 1\n", + "\\end{array}\n", + "\\end{array}\n", + "& \\begin{array}{c}\n", + "\\\\\n", + "\\begin{array}{c|cc}\n", + "0 & &\\\\\n", + "1/2 & 1/2 &\\\\ \\hline\n", + " & 0 & 1\n", + "\\end{array}\n", + "\\end{array}\n", + "& \\begin{array}{c|ccc}\n", + "0 & & & \\\\\n", + "1/3 & 1/3 & & \\\\\n", + "2/3 & 0 & 2/3 & \\\\ \\hline\n", + " & 1/4 & 0 & 3/4\n", + "\\end{array}\\\\\n", + "\\textrm{Euler} & \\textrm{Runge} & \\textrm{Heun}\n", + "\\end{array}$$\n", + "$$\\begin{array}[t]{cc}\n", + "\\begin{array}{c|cccc}\n", + "0 & & & &\\\\\n", + "1/2 & 1/2 & & &\\\\\n", + "1/2 & 0 & 1/2 & &\\\\ \n", + "1 & 0 & 0 & 1 &\\\\ \\hline\n", + " & 1/6 & 2/6 & 2/6 & 1/6\n", + "\\end{array}\n", + "& \\begin{array}{c|cc}\n", + "1/2-\\sqrt{3}/6 & 1/4 & 1/4-\\sqrt{3}/6\\\\\n", + "1/2+\\sqrt{3}/6 & 1/4+\\sqrt{3}/6 & 1/4\\\\ \\hline\n", + " & 1/2 & 1/2\n", + "\\end{array}\\\\\n", + "\\textrm{La méthode Rk4} & \\textrm{Gauss d'ordre 4}\n", + "\\end{array}$$\n" + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "metadata": {}, + "outputs": [], + "source": [ + "using LinearAlgebra\n", + "# import Pkg; Pkg.add(\"Printf\")\n", + "using Printf" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Les méthodes de Runge Kutta explicites" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "euler (generic function with 1 method)" + ] + }, + "execution_count": 19, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function euler(fun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real},\n", + " N::Integer)::Tuple{Vector{Float64},Array{Float64,2}}\n", + "# Numerical integration via Euler scheme\n", + "# Input\n", + "# -----\n", + "# fun : second member of the ode\n", + "# xpoint = fun(x,p,t)\n", + "# xpoint : value of fun\n", + "# array(n)\n", + "# x : state\n", + "# array(n)\n", + "# p : vector of parameter\n", + "# t : time\n", + "# real\n", + "# tspan : (t0, tf)\n", + "# x0 : initial condition\n", + "# array(n)\n", + "# p : vector of parameters\n", + "# N : number of step\n", + "# integer\n", + "# Output\n", + "# ------\n", + "# T : vector of time\n", + "# array(N+1)\n", + "# X : X[i,:] is the approximation of x(t_i)\n", + "# array(N+1,n)\n", + "# ----------------------------------------------\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "# A compléter\n", + " \n", + "return T,X\n", + "end\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "runge (generic function with 1 method)" + ] + }, + "execution_count": 20, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function runge(fun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real}, N::Integer)\n", + "# Numerical integration via the runge scheme\n", + "# the same as Euler\n", + "# ----------------------------------------------\n", + "\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "# A compléter\n", + " \n", + "return T,X\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "heun (generic function with 1 method)" + ] + }, + "execution_count": 21, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function heun(fun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real}, N::Integer)\n", + "# Numerical integration via the Heun scheme\n", + "# the same as Euler\n", + "# ----------------------------------------------\n", + "\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "# A compléter\n", + " \n", + "return T,X\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "rk4 (generic function with 1 method)" + ] + }, + "execution_count": 22, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function rk4(fun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real}, N::Integer)\n", + "# Numerical integration via the RK4 scheme\n", + "# the same as Euler\n", + "# ----------------------------------------------\n", + "\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "# A compléter\n", + " \n", + "return T,X\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Vérification des résultats dans le cas N = 10" + ] + }, + { + "cell_type": "code", + "execution_count": 23, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Euler\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + "Runge\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + "Heun\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + "RK4\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n", + " 0.000000000000000 0.000000000000000 0.000000000000000\n" + ] + } + ], + "source": [ + "function vdp(x::Vector{<:Real},p::Vector{<:Real},t::Real)\n", + "# Van der Pol model\n", + "# second member of the IVP\n", + "# p : paramter vector\n", + "# t : time, variable not use here\n", + "#\n", + "# A compéter\n", + " xpoint = similar(x)\n", + " return xpoint \n", + "end\n", + "\n", + "function affiche(T,X)\n", + "# Affiche les résulats : utile pour vérifier que l'on a les résultats attendus pour N=25\n", + "#\n", + " N = length(T)-1\n", + " for i in 1:N+1\n", + " @printf(\"%20.15f %20.15f %20.15f\\n\", T[i],X[i,1],X[i,2])\n", + " end\n", + "end\n", + "\n", + "mu = 1\n", + "p = [mu]\n", + "tf = 6.6632868593231301896996820305\n", + "tspan = (0.0, tf)\n", + "x0 = [2.00861986087484313650940188,0]\n", + "N = 10\n", + "\n", + "println(\"Euler\")\n", + "T_Euler,X_Euler = euler(vdp, tspan, x0, p, N)\n", + "affiche(T_Euler,X_Euler)\n", + "\n", + "println(\"Runge\")\n", + "T_Runge,X_Runge = runge(vdp, tspan, x0, p, N)\n", + "affiche(T_Runge,X_Runge)\n", + "\n", + "println(\"Heun\")\n", + "T_Heun,X_Heun = heun(vdp, tspan, x0, p, N)\n", + "affiche(T_Heun,X_Heun)\n", + "\n", + "println(\"RK4\")\n", + "T_RK4,X_RK4 = rk4(vdp, tspan, x0, p, N)\n", + "affiche(T_RK4,X_RK4)\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Graphiques de la figure 3, cas N = 25" + ] + }, + { + "cell_type": "code", + "execution_count": 24, + "metadata": {}, + "outputs": [], + "source": [ + "using Plots\n", + "#\n", + "# A compléter" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Ordre pour les méthodes explicites\n", + "On demande de réaliser la figure 1 uniquement pour la première variable" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": {}, + "outputs": [], + "source": [ + "# A compléter\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Runge Kutta implicite, schéma de Gauss\n", + "### Point fixe\n", + "#### Code\n" + ] + }, + { + "cell_type": "code", + "execution_count": 26, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "gauss_fp (generic function with 1 method)" + ] + }, + "execution_count": 26, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function gauss_fp(fun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real},\n", + " option)::Tuple{Vector{Float64},Array{Float64,2},Int,Array{Int,1}}\n", + "# Numerical integration with Gauss scheme and the fixed point algorithm\n", + "#\n", + "# Input parameters\n", + "# ----------------\n", + "# fun : second member of the ode\n", + "# fun(x,p,t)\n", + "# Input\n", + "# x : state\n", + "# array(n)\n", + "# p : vector of parameters\n", + "# t : time\n", + "# real\n", + "# Output\n", + "# xpoint : value of fun\n", + "# array(n) \n", + "# tspan = (t0,tf)\n", + "# x0 = initial point\n", + "# array(n)\n", + "# option[1] = N = number of step\n", + "# Integer\n", + "# option[2] = maximum number of iterations for the fixed point\n", + "# real\n", + "# option[3] = epsilon for the test of progress in the fixed point\n", + "# real\n", + "#\n", + "# Output parameters\n", + "# -----------------\n", + "# T = vector of time\n", + "# real(N+1)\n", + "# X = Matrix of solution\n", + "# real(N+1,n)\n", + "# The line i of [T X] contains ti and x_i the approximation of x(t_i)\n", + "# ifail[i] = number of iteration for the fixed point on [t_i,t_{i+1}] if the fixed point converge\n", + "# ifail[i] = -1 if the fixed point don't converge on [t_i,t_{i+1}]: maximum number of iteration\n", + "# is attained in the fixed point\n", + "# integer(N)\n", + "# nfun = number of evaluation of phi\n", + "# integer\n", + "# ndfun = number of evaluation of the derivative of fun\n", + "# integer\n", + "# ---------------------------------------------------------------------------------------------\n", + "#\n", + "# Initialisation\n", + "# --------------\n", + "# Coefficent of the Butcher's array\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "nfun = 0\n", + "ifail = ones(Int,N)\n", + "# A compléter\n", + "return T,X,nfun,ifail\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Test dans le cas N = 10, nb_itmax = 15 et K_eps = 1.e-12" + ] + }, + { + "cell_type": "code", + "execution_count": 27, + "metadata": {}, + "outputs": [], + "source": [ + "# A compléter" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Graphique 2" + ] + }, + { + "cell_type": "code", + "execution_count": 28, + "metadata": {}, + "outputs": [], + "source": [ + "# A compléter" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Ordre" + ] + }, + { + "cell_type": "code", + "execution_count": 29, + "metadata": {}, + "outputs": [], + "source": [ + "# A compléter" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Newton\n", + "#### Code" + ] + }, + { + "cell_type": "code", + "execution_count": 30, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "gauss_newton (generic function with 1 method)" + ] + }, + "execution_count": 30, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function gauss_newton(fun, dfun, tspan::Tuple{<:Real,<:Real}, x0::Vector{<:Real}, p::Vector{<:Real},\n", + " option)::Tuple{Vector{Float64},Array{Float64,2},Int,Int,Array{Int,1}}\n", + "# Numerical integration with Gauss scheme and the Newton algorithm\n", + "#\n", + "# Input parameters\n", + "# ----------------\n", + "# fun : second member of the ode\n", + "# fun(x,p,t)\n", + "# Input\n", + "# x : state\n", + "# array(n)\n", + "# p : vector of parameters\n", + "# t : time\n", + "# real\n", + "# Output\n", + "# xpoint : value of fun\n", + "# array(n) \n", + "# tspan = (t0,tf)\n", + "# x0 = initial point\n", + "# array(n)\n", + "# option[1] = N = number of step\n", + "# Integer\n", + "# option[2] = maximum number of iterations for the Newton algorithm\n", + "# real\n", + "# option[3] = epsilon for the test of progress in the Newton algorithm\n", + "# real\n", + "#\n", + "# Output parameters\n", + "# -----------------\n", + "# T = vector of time\n", + "# real(N+1)\n", + "# X = Matrix of solution\n", + "# real(N+1,n)\n", + "# The line i of [T X] contains ti and x_i the approximation of x(t_i)\n", + "# ifail[i] = number of iteration for the fixed point on [t_i,t_{i+1}] if the fixed point converge\n", + "# ifail[i] = -1 if the fixed point don't converge on [t_i,t_{i+1}]: maximum number of iteration\n", + "# is attained in the fixed point\n", + "# integer(N)\n", + "# nfun = number of evaluation of phi\n", + "# integer\n", + "# ndfun = number of evaluation of the derivative of fun\n", + "# integer\n", + "# ndfun = number of evaluation of the derivative of fun\n", + "# integer\n", + "# ---------------------------------------------------------------------------------------------\n", + "\n", + "n = length(x0) # dimension n\n", + "T = zeros(N+1)\n", + "X = zeros(N+1,n)\n", + "nfun = 0\n", + "ndfun = 0\n", + "ifail = ones(Int,N)\n", + "# A compléter\n", + "\n", + "return T,X,nfun,ndfun,ifail\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Test dans le cas N = 10, nb_itmax = 15 et K_eps = 1.e-12" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "dvdp (generic function with 1 method)" + ] + }, + "execution_count": 31, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "function dvdp(x::Vector{<:Real},p::Vector{<:Real},t::Real)::Matrix{<:Real}\n", + "# Van der Pol model\n", + "# derivative of the second member of the IVP\n", + "# x : state\n", + "# real(n)\n", + "# p : paramter vector\n", + "# t : time, variable not use here\n", + "# real\n", + "# Output\n", + "# Xpoint::Jacobian matrix{<:Real},\n", + "# real(n,n)\n", + "# --------------------------------------\n", + " n = length(x)\n", + " Jac = ones(n,n)\n", + "# A compléter\n", + " \n", + " return Jac\n", + "end\n", + "\n", + "# comptéter\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Graphique 2" + ] + }, + { + "cell_type": "code", + "execution_count": 32, + "metadata": {}, + "outputs": [], + "source": [ + "# A compléter\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Ordre" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": { + "scrolled": false + }, + "outputs": [], + "source": [ + "# A compléter\n" + ] + } + ], + "metadata": { + "celltoolbar": "Aucun(e)", + "kernelspec": { + "display_name": "Julia 1.4.1", + "language": "julia", + "name": "julia-1.4" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.4.1" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/ode/TP/sujet-tp.pdf b/ode/TP/sujet-tp.pdf new file mode 100644 index 0000000000000000000000000000000000000000..a07a6f398efcd0373d3756542d5ffca404fe28ce Binary files /dev/null and b/ode/TP/sujet-tp.pdf differ diff --git a/ode/cours_ode_ModIA/main.pdf b/ode/cours_ode_ModIA/main.pdf new file mode 100644 index 0000000000000000000000000000000000000000..dcf34157e6e76d359b8a6ea0020270bf55c83928 Binary files /dev/null and b/ode/cours_ode_ModIA/main.pdf differ