diff --git a/README.md b/README.md index ae73d50069c1118d4669a61d6c6f5912e0071cea..a4a82831a7714548cfd5010d3066e06abaa8eab7 100644 --- a/README.md +++ b/README.md @@ -63,6 +63,7 @@ Running `dpkg -i cassiopee_*.deb` will install Cassiopée in `/opt/Cassiopee` * mkdocs: `pip install mkdocs` * mkdocs-material: `pip install mkdocs-material` * python-markdown-math: `pip install https://github.com/mitya57/python-markdown-math/archive/master.zip` + * pygments: `pip install pygments` #### build .exe installer diff --git a/docs-fr/calculators/hsl/courbe_remous.md b/docs-fr/calculators/hsl/courbe_remous.md index 9b35577c3eacdef25d2672269a16a8dbaab0ceea..4ecaa51d277904b596f80f3831ea5731c1498656 100644 --- a/docs-fr/calculators/hsl/courbe_remous.md +++ b/docs-fr/calculators/hsl/courbe_remous.md @@ -1,37 +1,41 @@ -# Courbe de remous !! TODO !! +# Courbe de remous +Le calcul de la courbe de remous fait intervenir l'équation différentielle suivante : -Le régime uniforme se caractérise par une hauteur d'eau appelée hauteur normale. La hauteur normale est atteinte quand la ligne d'eau est parallèle au fond, la charge est alors elle-même parallèle à la ligne d'eau et donc la perte de charge est égale à la pente du fond : -\(I_f = J\) -Avec : +$$\frac{dy}{dx}=\frac{I_f - J(h)}{1-F^2(h)}$$ -- \(I_f\) : la pente du fond en m/m -- \(J\) : la perte de charge en m/m +où \(I_f\) est la pente d'un canal, \(J\) la formule nous donnant la perte de charge locale (dépendant la hauteur d'eau), \(y\) désigne ici la hauteur d'eau. -La perte de charge {J} est ici calculée avec la formule de Manning-Strickler : +On donne ainsi, pour un canal rectangulaire de largeur \(b\) et coefficient de Strickler \(K\): -$$J=\frac{U^2}{K^{2}R^{4/3}}=\frac{Q^2}{S^2K^{2}R^{4/3}}$$ +$$J=\frac{Q^2 (b+2y)^{4/3}}{K^2 b^{10/3}y^{10/3}}$$ -Avec : +et -- \(K\) : le coefficient de Strickler en m<sup>1/3</sup>/s +$$F^2=\frac{Q^2}{gb^2y^3}$$ -En régime uniforme, on obtient la formule : +L'intégration de l'équation pourra se faire par l'une des méthodes suivantes : [Runge-Kutta 4](../../methodes_numeriques/rk4.md), [Euler explicite](../../methodes_numeriques/euler_explicite.md), [intégration de trapèzes](../../methodes_numeriques/integration_trapezes.md). -$$Q=KR^{2/3}S\sqrt{I_f}$$ +En fonction du régime d'écoulement, le calcul pourra s'effectuer : -A partir de laquelle, on peut calculer analytiquement le débit \(Q\), la pente \(I_f\) et le Strickler \(K\) analytiquement. + * de l'aval vers l'amont pour le régime fluvial avec définition d'une condition limite aval + * de l'amont vers l'aval pour le régime torrentiel avec définition d'une condition limite amont -Pour calculer la hauteur normale \(h_n\) , on peut résoudre \(f(h_n)=Q-KR^{2/3}S\sqrt{I_f}=0\) +Si on prend l'exemple d'un canal rectangulaire, [l'exemple de code scilab proposé pour la résolution d'équation différentielle ordinaire](../../methodes_numeriques/euler_explicite.md) est modifié comme suit : -en utilisant la méthode de Newton : +```scilab + b=0.3; + K=50; + If=0.005; + Q=0.01; + function z=DQ(y);z=Q-K*(b*y)^(5/3)/(b+2*y)^(2/3)*sqrt(If); endfunction + yn=fsolve(0.5,DQ); + tmax=0; + t0=10; + dt=-0.5; + function z=f(y,t);z=(If-Q^2*(b+2*y)^(4/3)/(K^2*(b*y)^(10/3)))/(1-Q^2/(9.81*b^2*y^3)); endfunction + y0=0.12; +``` -$$h_{k+1} = h_k - \frac{f(h_k)}{f'(h_k)}$$ - - avec : - -- \(f(h_k) = Q-KR^{2/3}S\sqrt{I_f}\) -- \(f'(h_k) = -K \sqrt{I_f}(\frac{2}{3}R'R^{-1/3}S+R^{2/3}S')\) - -Pour calculer les paramètres géométriques de la section, le module de calcul utilise l'équation de calcul du débit et résout le problème par dichotomie. +ce qui nous donne la profondeur normale, et le tracé de la ligne d'eau. En fonction de la méthode numérique utilisée, on peut avoir de grosses erreurs dans le cas d'une courbe de remous de type F2 (condition aval sous la hauteur normale), car les pentes de ligne d'eau sont beaucoup plus fortes, et donc beaucoup plus sujettes aux erreurs liées à l'interpolation linéaire. On en déduit donc que d'une part le choix de la méthode de résolution est important, d'autre part qu'il est essentiel d'avoir un regard critique sur les solutions (avec une interprétation des processus qu'on cherche à modéliser). diff --git a/docs-fr/methodes_numeriques/brent.md b/docs-fr/methodes_numeriques/brent.md new file mode 100644 index 0000000000000000000000000000000000000000..0b72a67a240bc6a5e6c5d565a6579ec2f2f039c6 --- /dev/null +++ b/docs-fr/methodes_numeriques/brent.md @@ -0,0 +1 @@ +# Méthode de Brent diff --git a/docs-fr/methodes_numeriques/euler_explicite.md b/docs-fr/methodes_numeriques/euler_explicite.md new file mode 100644 index 0000000000000000000000000000000000000000..583035de39167f252b8a5623d7d8144be25377f3 --- /dev/null +++ b/docs-fr/methodes_numeriques/euler_explicite.md @@ -0,0 +1,50 @@ +# Méthode d'Euler explicite + +Pour décrire un processus d'évolution, ou le profil d'une ligne d'eau, par exemple, on est souvent amené à résoudre une équation différentielle ordinaire (EDO) du premier ordre. Cette équation écrit comment varie une fonction, en un point donné (un instant ou un point de l'espace), connaissant la valeur de cette fonction mathématique, le problème à résoudre s'écrit: + +$$ +\left\{ + \begin{array}{rcr} + \frac{dy}{dt} & = & f(y,t) \\ + y(t=t_0) & = & y_0 \\ + \end{array} +\right. +$$ + +où \(\frac{dy}{dt}\) désigne la dérivée par rapport à t de la fonction \(y\) (qui dépend de la variable \(t\)); la variable \(y_0\) est appelée la condition à la limite; elle conditionne la solution finale de l'équation. + +Comme souvent on ne connait pas de solution analytique de ce problème, on va utiliser des méthodes approchées pour estimer la solution. On fait donc une discrétisation de la variable \(t\). On note ainsi \(\Delta t\) le pas de discrétisation, et on résout le problème aux points \(t_0\), \(t_1=t_0+\Delta t\), \(t_2=t_0+2\Delta t\), ..., \(t_n=t_0+n\Delta t\) où \(n\) est un entier. + + +La méthode d'Euler explicite est la plus intuitive; elle consiste à considérer que, d'un point \(t_i\) au point \(t_{i+1}\), la fonction évolue linéairement, avec une trajectoire qui est celle qu'on peut calculer au point \(t_i\). + +Le problème se résout donc de la façon suivante: + + * on connait la fonction \(f\), un point \(t_i\) où on connait \(y_i\) + * on peut donc calculer \(y'_i=f(y,t)\) + * on estime alors la valeur de \(y\) au point \(t_{i+1}=t_i+\Delta t\) : \(y_{i+1}\simeq y_i + y'_i \Delta t\) + * on peut alors itérer (résoudre pas à pas) pour passer au point suivant. Le problème est initialisé en partant de \(t_0\) où on connait \(y_0\) (condition à la limite). + +On sent bien que ce schéma pourra donner de bons résultats uniquement si \(\Delta t\) n'est pas trop grand. Des valeurs de \(\Delta t\) trop grandes peuvent donner des résultats complètement faux, conduisant à des interprétations physiques erronées. Son intérêt est toutefois sa simplicité, et il s'implémente facilement sur un tableau. + + +## Exemple d'application: processus exponentiel + +Considérons le problème (simple) suivant: +$$ +\left\{ + \begin{array}{rcr}\label{eq:exp} + \frac{dy}{dt} & = & -a y \\ + y(t=t_0) & = & y_0 \\ + \end{array} +\right. +$$ + +On a donc ici \(f(y,t)=-ay\). La solution analytique se résout facilement, donnant \(y(t)=y_0 \exp\left(-a(t-t_0)\right)\). +On peut résoudre le problème par la méthode d'Euler: + + * on choisit \(\Delta t\) (par exemple, \(\Delta t=1\)) + * calculer \( y_1=y_0 - a y_0 \Delta t\) + * calculer \( y_2=y_1 - a y_1 \Delta t\) etc. + +On constate que la résolution n'est pas très précise; ceci est lié au pas de calcul trop grand compte tenu de la méthode choisie et de l'équation à résoudre. diff --git a/docs-fr/methodes_numeriques/integration_trapezes.md b/docs-fr/methodes_numeriques/integration_trapezes.md new file mode 100644 index 0000000000000000000000000000000000000000..7616f96eb7e4d35750d833fd04f3cfd1afe44fdb --- /dev/null +++ b/docs-fr/methodes_numeriques/integration_trapezes.md @@ -0,0 +1 @@ +# Méthode par intégration de trapèzes diff --git a/docs-fr/methodes_numeriques/newton.md b/docs-fr/methodes_numeriques/newton.md new file mode 100644 index 0000000000000000000000000000000000000000..040e4c8e044033d449900f663e82f4ab111581e9 --- /dev/null +++ b/docs-fr/methodes_numeriques/newton.md @@ -0,0 +1 @@ +# Méthode de Newton diff --git a/docs-fr/methodes_numeriques/rk4.md b/docs-fr/methodes_numeriques/rk4.md new file mode 100644 index 0000000000000000000000000000000000000000..202c488ce7ccb0d542075d1891622de4227b71f3 --- /dev/null +++ b/docs-fr/methodes_numeriques/rk4.md @@ -0,0 +1,64 @@ +# Schéma de Runge-Kutta d'ordre 4 + +Le schéma de Runge-Kutta d'ordre 4 est basée sur une approximation de la dérivée à un ordre supérieur (ordre 4). Le principe de la discrétisation est le même que pour la [méthode d'Euler](euler_explicite.md), mais on va faire quelques calculs supplémentaires pour approcher la dérivée: + + * on connait la fonction \(f\), un point \(t_i\) où on connait \(y_i\) + * on peut donc calculer \(f_1=y'_i=f(y_i,t_i)\) (cf. méthode d'Euler=pente au point (\(t_i,y_i\))) + * on calcule \(f_2=f(y_i+\frac12 \Delta t f_1,t_i+\frac12 \Delta t)\) (valeur estimée au milieu de l'intervalle, avec la pente prise en \(t_i\)) + * on calcule \(f_3=f(y_i+\frac12 \Delta t f_2,t_i+\frac12 \Delta t)\) (valeur estimée au milieu de l'intervalle \(t_{i+1/2}\), avec la pente prise en \(t_{i+1/2}\)) + * on calcule \(f_4=f(y_i+ \Delta t f_3,t_i+ \Delta t)\) (valeur estimée en \(t_{i+1}\), avec la pente prise en \(t_{i+1/2}\)) + * on a alors \(y_{i+1}\simeq y_1 + \frac16 \Delta t(f_1+2f_2+2f_3+f_4)\) + * on peut alors itérer (résoudre pas à pas) pour passer au point suivant. Le problème est initialisé en partant de \(t_0\) où on connait \(y_0\) (condition à la limite). + +On voit clairement que la méthode est beaucoup plus précise. Même avec un pas de calcul beaucoup plus élevé, on approche correctement la solution, alors que la méthode d'Euler donne des résultats très éloignés de la solution exacte. On remarque aussi que, entre les points de discrétisation, on a approché la solution par des segments de droite (bien qu'on puisse avoir une interpolation plus évoluée si on a besoin de connaître des valeurs intermédiaires). + +Le programme peut s'écrire de la façon suivante, avec le logiciel de calcul scientifique Scilab: + +```scilab +// Programme de resolution d'une equation differentielle +// Equation a resoudre: dy/dt=f(y,t) + +// Definition de la fonction f +a=-0.1; +function z=f(y,t);z=a*y; endfunction + +// Condition à la limite: +t0=0; +y0=1; + +// Discretisation du temps +tmax=50; +dt=5; + +// Nombre de pas de discretisation +N=(tmax-t0)/dt; +// Indices +ii=1:N+1; +t=(ii-1)*dt; // vecteur temps + +t2=0:tmax; // vecteur temps avec un pas plus fin, pour la solution exacte + +// Solution par méthode d'Euler, notée ye +ye(1)=y0; // condition à la limite +for i=1:N + ye(i+1)=ye(i)+dt*f(ye(i),t(i)); +end + +// Solution par méthode RK4, notée yrk4 +yrk4(1)=y0; // condition à la limite +for i=1:N + f1=f(yrk4(i),t(i)); + f2=f(yrk4(i)+f1*dt/2,t(i)+dt/2); + f3=f(yrk4(i)+f2*dt/2,t(i)+dt/2); + f4=f(yrk4(i)+f3*dt,t(i)+dt); + yrk4(i+1)=yrk4(i)+dt*(f1+2*f2+2*f3+f4)/6; +end + +// Tracé des solutions +scf(2) +plot(t2,exp(a*t2),'k-',t,ye,'kd-',t,yrk4,'ks-') +title('Resolution par méthode de Runge-Kutta d''ordre 4 - dy/dt=-0.1y') +xlabel('t') +ylabel('y') +legend('Solution exacte','Méthode Euler, dt=5','Méthode RK4, dt=5') +``` \ No newline at end of file diff --git a/mkdocs.yml b/mkdocs.yml index 1353c78018d3e1604c0209f0712f6f7ed9cf7946..c8a9621dd2bdb6d04e05f5cedb47f9327c632c75 100644 --- a/mkdocs.yml +++ b/mkdocs.yml @@ -14,6 +14,7 @@ extra_javascript: markdown_extensions: - mdx_math - footnotes + - codehilite nav: - avant-propos.md - Présentation de Cassiopée: @@ -25,7 +26,7 @@ nav: - calculators/hyd_en_charge/lechapt-calmon.md - calculators/hyd_en_charge/cond_distri.md - Hydraulique à surface libre: - - calculators/hsl/regime_uniforme.md + - Régime uniforme: calculators/hsl/regime_uniforme.md - calculators/hsl/courbe_remous.md - calculators/hsl/section_parametree.md - calculators/hsl/types_sections.md @@ -54,5 +55,11 @@ nav: - Dévalaison: - Perte de charge sur grille de prise d'eau: calculators/devalaison/grille.md - Impact de jet: calculators/devalaison/jet.md + - Méthodes numériques de résolution : + - Runge-Kutta 4: methodes_numeriques/rk4.md + - Euler explicite: methodes_numeriques/euler_explicite.md + - Intégration de trapèzes: methodes_numeriques/integration_trapezes.md + - Méthode de Brent: methodes_numeriques/brent.md + - Méthode de Newton: methodes_numeriques/newton.md - CHANGELOG.md - mentions_legales.md