...

Commits (3)
 ... ... @@ -110,14 +110,14 @@ \end{equation*} De plus, à l'aide des méthodes de {\it Jacobi} et {\it Gauss-Seidel} (voir Série 2) essayez d'obtenir la figure suivante. \begin{figure}[H] \centering \includegraphics[scale=0.6]{convergence.pdf} \caption{La ligne bleue correspond à la méthode de Jacobi, la ligne rouge correspond à la méthode de Gauss-Seidel et la ligne violette correspond à la méthode SOR. Les lignes noires pointillées représentent les estimations de la convergence des trois méthodes.} \end{figure} % \begin{figure}[H] % \centering % \includegraphics[scale=0.6]{convergence.pdf} % \caption{La ligne bleue correspond à la méthode de Jacobi, %la ligne rouge correspond à la méthode de Gauss-Seidel %et la ligne violette correspond à la méthode SOR. %Les lignes noires pointillées représentent les estimations de la convergence des trois méthodes.} % \end{figure} \end{enumerate} ... ...
 AAltermatt Anne-CatherineBerisV CanosaAlvesL CorviniV FolL GeimerA GrosjeanL HorvathR MaizCalleP Perrin-GanierN PistoresiL RivollatQ RohrbasserL StraubC WahlN
 ... ... @@ -5,22 +5,32 @@ function A=FDLaplacian(n,d) h=1/(n+1); D=1/h^2*spdiags(ones(n,1)*[1 -2*d 1],[-1 0 1],n,n); D=1/h^2*spdiags(ones(n,1)*[1 -2 1],[-1 0 1],n,n); if d==1 A=D; A=D; elseif d==2 X=spdiags(ones(n,1),0,n,n); I2=1/h^2*spdiags(ones(n^2,1)*[1 1],[-n n],n^2,n^2); A=kron(X,D)+I2; I = eye(n); A=kron(I,D) + kron(D,I); elseif d==3 X=spdiags(ones(n^2,1),0,n^2,n^2); I2=1/h^2*spdiags(ones(n^3,1)*[1 1],[-n n],n^3,n^3); I3=1/h^2*spdiags(ones(n^3,1)*[1 1],[-n^2 n^2],n^3,n^3); A=kron(X,D)+I2+I3; I = eye(n); D1=kron(I,D); D2=kron(D,I); II=kron(I,I); A=kron(D1,I) + kron(D2,I) + kron(II,D); end % if d==1 % A=D; % elseif d==2 % X=spdiags(ones(n,1),0,n,n); % I2=1/h^2*spdiags(ones(n^2,1)*[1 1],[-n n],n^2,n^2); % A=kron(X,D)+I2; % elseif d==3 % X=spdiags(ones(n^2,1),0,n^2,n^2); % I2=1/h^2*spdiags(ones(n^3,1)*[1 1],[-n n],n^3,n^3); % I3=1/h^2*spdiags(ones(n^3,1)*[1 1],[-n^2 n^2],n^3,n^3); % A=kron(X,D)+I2+I3; % end % avec blkdiag ... ...
 % Series 1 Test Suite n=101; d=1; T=FDLaplacianTrue(n-2,d); A=FDLaplacien(n-2,d); disp(normest(A-T,2)); disp(normest((n-1)^2 * A-T,2)); % x=linspace(0,1,n)'; f=sin(x(2:end-1)); % error=abs(A*f + f); % disp(norm(error(2:end-1),2)) n=11; d=2; T=FDLaplacianTrue(n-2,d); A=FDLaplacien(n-2,d); disp(normest(A-T,2)); disp(normest((n-1)^2 * A - T,2)); % x=linspace(0,1,n); % y=linspace(0,1,n)'; % f=sin(y(2:end-1)) * sin(x(2:end-1)); % error=abs(A*f(:) + f(:)); % disp(norm(error(2:end-1),2)); n=7; d=3; T=FDLaplacianTrue(n-2,d); A=FDLaplacien(n-2,d); disp(normest(A-T,2)); disp(normest((n-1)^2 * A - T,2)); \ No newline at end of file
 function u=Heat2D(f,g,p) function [u,F]=Heat2D(f,g,p) % HEAT2D solves the stationary heat equation over a square % u=Heat2D(f,g,p) solves the heat equation over a square of size one % by one. The source function is given by f, the boundary condition is % given by g, the number of grid points for each direction is a quarter % of the length(g) and one chooses the iterative solver with p, % i.e. p=0 uses Jacobi, p=1 uses Gauss-Seidel. n = length(g); n = n/4; ind=1:n; I = eye(n); G = kron(g(ind),I(:,1)) + kron(g(ind+n),I(:,end)) + ... kron(I(:,1),g(ind+2*n)) + kron(I(:,end),g(ind+3*n)); G = (n+1)^2 * G; F = f - G; if p==0 LU = (n+1)^2*spdiags(ones(n,1)*[1 0 1],[-1 0 1],n,n); LU = kron(I,LU) + kron(LU,I); D = 1/( 4 * (n+1)^2 ); elseif p==1 end u0 = kron(zeros(n,1),zeros(n,1)); err= 1; tol= 1e-8; itermax = 1000; iter = 1; while err > tol && iter < itermax u = D*LU*u0 - D*F; err = norm(u-u0); u0= u; iter = iter+1; end \ No newline at end of file
 function D=Laplacian1D(n) % LAPLACIAN1D compute a finite difference Laplacian % A=FDLaplacian(n,d) computes a finite difference Laplacian on the % unit interval using n interior points: 0 .. h .. 2h .. .. nh .. 1 h=1/(n+1); D=1/h^2*spdiags(ones(n,1)*[1 -2 1],[-1 0 1],n,n); \ No newline at end of file
 % test suite for serie 2 n = 10; g = zeros(4*n,1); f = -ones(n^2,1); [u,F] = Heat2D(f,g,0); u = reshape(u,n,n); surf(u) L = FDLaplacianTrue(n,2); L*u(:) - F \ No newline at end of file
 ... ... @@ -10,6 +10,7 @@ \addtolength{\voffset}{-2cm} \usepackage{hyperref} \usepackage{float} %% commandes utiles %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \newcommand{\R}{\mathbb{R}} ... ... @@ -33,7 +34,7 @@ \begin{enumerate} \item (2 pts) Ecrire les méthodes itératives de {\it Jacobi} et {\it Gauss-Seidel} sous les deux formes standard, de correction, de résidu et de différences. sous les formes standard, de correction, de résidu et de différences. \item (4 pts) Soit $\vec{B}\in\C^{n\times n}$ une matrice carrée dont le rayon spectral $\rho(\vec{B})$ vérifie $\rho(\vec{B})<1$. Démontrer que ... ... @@ -74,7 +75,7 @@ En Matlab, la fonction \texttt{kron(A,B)} calcule le produit de Kronecker. \begin{enumerate} \item Construire le Laplacien discret en dimension deux en utilisant \texttt{kron}. Montrer qu'il est le m\^eme du r\'esultat de la fonction \texttt{FDLaplacian(n,2)} d'exercice 5 de s\'erie 1. Montrer qu'il est le m\^eme du r\'esultat de la fonction \texttt{FDLaplacian(n,2)} de l'exercice 5 de la s\'erie 1. \item Soit $\left \{ \lambda_i \right \}$ l'ensemble des valeurs propres du Laplacien discret en dimension une. D\'emontrer que $\left \{ \lambda_i + \lambda_j \right \}$ est l'ensemble des valeurs propres du Laplacien discret en dimension deux. Quelle est la relation entre les vecteurs propres? ... ...
 \documentclass[a4paper]{article} \input{style/head.tex} %------------------------------- % TITLE VARIABLES (identify your work!) %------------------------------- \newcommand{\yourname}{} % replace YOURNAME with your name \newcommand{\yournetid}{} % replace YOURNETID with your NetID \newcommand{\youremail}{} % replace YOUREMAIL with your email \newcommand{\assignmentnumber}{2} % replace X with assignment number \usepackage{preamble} \usepackage{pdfpages} \theoremstyle{definition} \newtheorem{ex}{Exercice} \begin{document} %------------------------------- % TITLE SECTION (do not modify unless you really need to) %------------------------------- \input{style/header.tex} %------------------------------- % ASSIGNMENT CONTENT (add your responses) %------------------------------- \begin{ex} Puisque $\rho(B) < 1$, $1 \notin \sigma(B)$, o\u $\sigma(B)$ est le spectre de $B$. Alors, $\det(I-B)\neq 0$, qui veut dire que $(I-B)$ est inversible. To show the Neumann series, we first consider the partial sum $\sum_{j=0}^n B^j$ and multiply it by $(I-B)$ to get \begin{equation*} (I-B) \sum_{j=0}^n B^j = \sum_{j=0}^n B^j - B^{j+1} = B^0 - B^{n+1} = I - B^{n+1}, \end{equation*} where we used the fact that the second term is a telescoping sum. Now, we multiply both terms by $(I-B)^{-1}$ to obtain \begin{equation*} \sum_{j=0}^n B^j = (I-B)^{-1} ( I - B^{n+1} ) = (I-B)^{-1} - (I-B)^{-1} B^{n+1}. \end{equation*} By subtracting $(I-B)^{-1}$, we get \begin{equation*} \sum_{j=0}^n B^j - (I-B)^{-1} = - (I-B)^{-1} B^{n+1}. \end{equation*} For any arbitrary $\epsilon>0$ and $k$ large enough we have that $\| B^k \| \leq ( \rho(B) + \epsilon )^k$. Now, we write \begin{equation*} \| \sum_{j=0}^n B^j - (I-B)^{-1} \| = \| (I-B)^{-1} B^{n+1} \| \leq \| (I-B)^{-1} \| \| B^{n+1} \| \leq \| (I-B)^{-1} \| ( \rho(B) + \epsilon )^{n+1}. \end{equation*} Since $\rho(B)<1$, we can choose $\epsilon$ sufficiently small such that $( \rho(B) + \epsilon ) < 1$ holds.\footnote{Notice also that $\rho(B)<1$, we can find a matrix norm $\| \cdot \|$ such that $\| B \| < 1$ (see Problem \ref{prob:2.3}), therefore we have $\|B^k\|\leq\|B\|^k$ and $\lim_{k\rightarrow\infty}\|B\|^k = 0$.} Taking the limit for $n \rightarrow \infty$ we get \begin{equation*} \lim_{n \rightarrow \infty} \| \sum_{j=0}^n B^j - (I-B)^{-1} \| \leq \lim_{n \rightarrow \infty} \| (I-B)^{-1} \| ( \rho(B) + \epsilon )^{n+1} = 0. \end{equation*} This means that \begin{equation*} \lim_{n \rightarrow \infty} \sum_{j=0}^n B^j = (I-B)^{-1}. \end{equation*} \end{ex} \begin{ex} J: M=D, N=A-D. GS: M=D+L, N=-U. For a splitting $A = M - N$ we have that \begin{eqnarray*} \vec{u}_{k+1} &=& M^{-1}N \vec{u}_k + M^{-1}\vec{f} \qquad \quad \; \; \: \text{(standard form)}, \nonumber \\ \vec{u}_{k+1} &=& \vec{u}_k + M^{-1}\vec{r}_k \qquad \qquad \qquad \; \, \text{(correction form)}, \nonumber \\ \vec{e}_{k+1} &=& M^{-1}N \vec{e}_k \qquad \qquad \qquad \qquad \text{(error form)}, \nonumber \\ \vec{d}_{k+1} &=& M^{-1}N \vec{d}_{k} \qquad \qquad \qquad \qquad \text{(difference form)}, \nonumber \\ \vec{r}_{k+1} &=& (I-AM^{-1})\vec{r}_{k} = N M^{-1} \vec{r}_{k} \; \text{(residual form)}. \nonumber \\ \end{eqnarray*} Now, it is sufficient to replace into these formulas the splitting matrices $M$ and $N$ corresponding to the methods of Jacobi, Gauss-Seidel, and SOR. For the method of Jacobi, we have \begin{eqnarray*} \vec{u}_{k+1} &=& -D^{-1}(L+U) \vec{u}_k + D^{-1}\vec{f} , \nonumber \\ \vec{u}_{k+1} &=& \vec{u}_k + D^{-1}\vec{r}_k , \nonumber \\ \vec{e}_{k+1} &=& -D^{-1}(L+U) \vec{e}_k , \nonumber \\ \vec{d}_{k+1} &=& -D^{-1}(L+U) \vec{d}_{k} , \nonumber \\ \vec{r}_{k+1} &=& (I-AD^{-1})\vec{r}_{k} = -(L+U)D^{-1} \vec{r}_{k} . \nonumber \\ \end{eqnarray*} For the method of Gauss-Seidel, we have \begin{eqnarray*} \vec{u}_{k+1} &=& -(D+L)^{-1}U \vec{u}_k + (D+L)^{-1}\vec{f} , \nonumber \\ \vec{u}_{k+1} &=& \vec{u}_k + (D+L)^{-1}\vec{r}_k , \nonumber \\ \vec{e}_{k+1} &=& -(D+L)^{-1}U \vec{e}_k , \nonumber \\ \vec{d}_{k+1} &=& -(D+L)^{-1}U \vec{d}_{k} , \nonumber \\ \vec{r}_{k+1} &=& (I-A(D+L)^{-1})\vec{r}_{k} = -U(D+L)^{-1} \vec{r}_{k} . \nonumber \\ \end{eqnarray*} Finally, for the SOR method, we have \begin{eqnarray*} \vec{u}_{k+1} &=& (\frac{1}{\omega}D+L)^{-1}(-U+(\frac{1}{\omega}-1)D) \vec{u}_k + (\frac{1}{\omega}D+L)^{-1}\vec{f}, \nonumber \\ \vec{u}_{k+1} &=& \vec{u}_k + (\frac{1}{\omega}D+L)^{-1}\vec{r}_k , \nonumber \\ \vec{e}_{k+1} &=& (\frac{1}{\omega}D+L)^{-1}(-U+(\frac{1}{\omega}-1)D) \vec{e}_k , \nonumber \\ \vec{d}_{k+1} &=& (\frac{1}{\omega}D+L)^{-1}(-U+(\frac{1}{\omega}-1)D) \vec{d}_{k} , \nonumber \\ \vec{r}_{k+1} &=& (I-A(\frac{1}{\omega}D+L)^{-1})\vec{r}_{k} = (-U+(\frac{1}{\omega}-1)D)(\frac{1}{\omega}D+L)^{-1} \vec{r}_{k} , \nonumber \\ \end{eqnarray*} that are equivalent to \begin{eqnarray*} \vec{u}_{k+1} &=& (D+\omega L)^{-1}(-\omega U+(1-\omega)D) \vec{u}_k + \omega(D+\omega L)^{-1}\vec{f}, \nonumber \\ \vec{u}_{k+1} &=& \vec{u}_k + \omega (D+\omega L)^{-1}\vec{r}_k , \nonumber \\ \vec{e}_{k+1} &=& (D+\omega L)^{-1}(-\omega U+(1-\omega)D) \vec{e}_k , \nonumber \\ \vec{d}_{k+1} &=& (D+\omega L)^{-1}(-\omega U+(1-\omega)D) \vec{d}_{k} , \nonumber \\ \vec{r}_{k+1} &=& (I-\omega A(D+\omega L)^{-1})\vec{r}_{k} = (-\omega U+(1-\omega)D)(D+\omega L)^{-1} \vec{r}_{k} . \nonumber \\ \end{eqnarray*} \end{ex} \begin{ex} La convergence arrive si $\norm{M^{-1}N}<1$. Toutes les normes de matrices sont \'equivalentes. \begin{align*} \norm{M^{-1}N}_\infty & = \max_i \sum_{j \neq i} \frac{\abs{a_{ij}}}{a_{ii}} < 1 \end{align*} \end{ex} \newcommand{\G}{\begin{bmatrix} 5/2 & -1 \\ 1 & 0 \end{bmatrix}} \begin{ex} \begin{align*} v_1 & = \G \begin{bmatrix} 1 \\ 2 \end{bmatrix} \\ & = \begin{bmatrix} 1/2 \\ 1 \end{bmatrix} \\ v_2 & = \G \begin{bmatrix} 1/2 \\ 1 \end{bmatrix} \\ & = \begin{bmatrix} 1/4 \\ 1/2 \end{bmatrix} \\ \det \left (\G - \lambda I \right ) & = \lambda^2 - \frac{5}{2}\lambda + 1 = 0 \\ \implies \lambda & = \frac{5/2 \pm \sqrt{25/4 - 16/4}}{2} \\ & = \frac{5/2 \pm 3/2}{2} \\ & = \begin{cases} 2 \\ 1/2 \end{cases} \\ \left (\G - 2I \right ) v = 0 & \implies v = \alpha \begin{bmatrix} 2 \\ 1 \end{bmatrix} \\ \left (\G - 1/2 I \right ) v = 0 & \implies v = \alpha \begin{bmatrix} 1/2 \\ 1 \end{bmatrix} \end{align*} \end{ex} %------------------------------------------------ \end{document}
 % Master Preamble \NeedsTeXFormat{LaTeX2e} \ProvidesPackage{preamble} \RequirePackage{amsmath} \RequirePackage{amssymb} \RequirePackage{amsthm} \RequirePackage{graphicx} \RequirePackage{hyperref} \RequirePackage{subcaption} \newcommand{\dxdy}{\frac{d #1}{d #2}} \newcommand{\dxdyk}{\frac{d^{#3} #1}{d {#2}^{#3}}} \newcommand{\pdxdy}{\frac{\partial #1}{\partial #2}} \newcommand{\liminfty}{\lim_{#1 \to \infty}} \newcommand{\limab}{\lim_{#1 \to #2}} \newcommand{\abs}{\left \vert #1 \right \vert} \newcommand{\norm}{\left \Vert #1 \right \Vert} \newcommand{\order}{\mathcal{O} \left ( #1 \right )} \newcommand{\set}{\left \{ #1 \right \}} \newcommand{\Set}{\left \{ #1 \ \middle \vert \ #2 \right \}} \newcommand{\vmat}{\begin{vmatrix} #1 \end{vmatrix}} \DeclareMathOperator{\sign}{sign} \newcommand{\bbn}{\mathbb{N}} \newcommand{\bbz}{\mathbb{Z}} \newcommand{\bbq}{\mathbb{Q}} \newcommand{\bbr}{\mathbb{R}} \newcommand{\bbc}{\mathbb{C}} \newcommand{\bbf}{\mathbb{F}} \newtheorem{lemma}{Lemma} \newtheorem{thm}{Theorem} \newtheorem{cor}{Corollary} \newtheorem{prop}{Proposition} \newtheorem{conj}{Conjecture} \newtheorem{defn}{Definition} \endinput \ No newline at end of file
 \fancyhead[C]{} \hrule \medskip \begin{minipage}{0.295\textwidth} \raggedright \footnotesize \yourname \hfill\\ \yournetid \hfill\\ \youremail \end{minipage} \begin{minipage}{0.4\textwidth} \centering \large S\'{e}rie \assignmentnumber\\ \end{minipage} \begin{minipage}{0.295\textwidth} \raggedleft \hfill\\ \end{minipage} \medskip\hrule \bigskip
 % Source: ss17_wissschreib (Eva) \lstset{ basicstyle=\ttfamily\scriptsize\mdseries, keywordstyle=\bfseries\color[rgb]{0.171875, 0.242188, 0.3125}, identifierstyle=, commentstyle=\color[rgb]{0.257813, 0.15625, 0}, stringstyle=\itshape\color[rgb]{0.0195313, 0.195313, 0.0117188}, numbers=left, numberstyle=\tiny, stepnumber=1, breaklines=true, frame=none, showstringspaces=false, tabsize=4, backgroundcolor=\color[rgb]{0.98,0.98,0.98}, captionpos=b, float=htbp, language=Python, xleftmargin=15pt, xrightmargin=15pt } %(deutsche) Sonderzeichen \lstset{literate=% {Ä}{{\"A}}1 {Ö}{{\"O}}1 {Ü}{{\"U}}1 {ä}{{\"a}}1 {ö}{{\"o}}1 {ü}{{\"u}}1 {ß}{{\ss}}1 } %Verwendung im Text: %-> \begin{lstlisting}[language=Python,firstnumber=27] ... \end{lstlisting} %-> \begin{lstlisting}[language=Python,numbers=none] ... \end{lstlisting} %-> \lstinline[language=JAVA]{...} \ No newline at end of file
 \documentclass[12pt,a4paper]{article} \usepackage{amsmath,amssymb,amsthm,epsfig,amscd,verbatim} \usepackage[french]{babel} \usepackage{fancyvrb} %\usepackage[T1]{fontenc} \usepackage[utf8]{inputenc} \addtolength{\textwidth}{3cm} \addtolength{\textheight}{4cm} \addtolength{\hoffset}{-2cm} \addtolength{\voffset}{-2cm} \usepackage{hyperref} \usepackage{float} %% commandes utiles %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \newcommand{\R}{\mathbb{R}} \newcommand{\C}{\mathbb{C}} \renewcommand{\vec}{\mathbf{#1}} \newcommand{\tr}{\mbox{tr}} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{document} \noindent UNIVERSIT\'E DE GEN\EVE \hfill Section de Math\'ematiques \\ \noindent Facult\'e des Sciences \hfill \A rendre le mercredi 1 avril \\[-3mm] \hrule \bigskip \begin{center} \textbf{\Large{Méthodes Itératives}} \\[1mm]\bigskip \textbf{Série 3} \end{center} \bigskip \begin{enumerate} \item \textbf{M\'ethodes it\'eratives stationnaires} \begin{enumerate} \item (3 pts) Considérer le Laplacien discret $A$ de $-\Delta$ en dimension 2. Démontrer que la méthode de Jacobi et la méthode de Richardson (avec le paramètre optimal $\alpha_{opt}$) ont le même comportement de convergence. \textit{Indice:} Pour trouver $\alpha_{opt}$, rappeler exercice 5(b) de la s\'erie 2 et que les valeurs propres d'une matrice tridiagonale et Toeplitz de taille $n \times n$ sont \begin{equation*} \lambda_n = a - 2 \sqrt{bc} \cos \left ( \frac{k \pi}{n+1} \right ), \ k = 1,...,n \end{equation*} o\u $a$, $b$ et $c$ sont les valeurs sur les diagonales (primaire, supérieure et inférieure, respectivement). \item (3 pts) Considérer une matrice $A \in \R^{2 \times 2}$ et le splitting $A = D + L + U$ avec $D$ inversible. Démontrer que si la méthode de Jacobi converge, alors la méthode de Gauss-Seidel converge aussi et plus rapidement. \item (4 pts) Considérer le Laplacien discret $A\in \R^{n^2\times n^2}$ de $-\Delta$ en dimension 2 et soit $h=\frac{1}{n+1}$ la taille de la grille. Dénoter avec $\rho(G_{\rm J};h)$ et $\rho(G_{\rm GS};h)$ les rayons spectraux des méthodes de Jacobi et Gauss-Seidel en fonction de $h$. Démontrer que : \begin{itemize}%\itemsep0em \item $\rho(G_{\rm J};h) = 1 - \frac{\pi^2}{2} h^2 + O(h^4)$. \item $\rho(G_{\rm GS};h) = 1 - \pi^2 h^2 + O(h^4)$. \end{itemize} \end{enumerate} \item \textbf{M\'ethode de descente maximale} \begin{enumerate} \item (2 pts) Consid\'erer la fonction $F: \R^n \rightarrow \R$ d\'efinit par $F(u) = \frac{1}{2} u^\top A u - u^\top f$ pour une $A \in \R^{n \times n}$ sym\'etrique et $f \in \R^n$. D\'emontrer que $\nabla F(u) = Au - f$. \item (4 pts) Considérer une fonction $f : \R^n \rightarrow \R$ et on dénote avec $\nabla f(\vec{v})$ le gradient de $f$ en $\vec{v} \in \R^n$ obtenu au moyen du produit scalaire habituel pour $\R^n$. Démontrer que $\nabla f(\vec{v})$ est orthogonal à la ligne de niveau de $f$ en $\vec{v}$ et que $-\nabla f(\vec{v})$ est la direction de descente maximale de $f$ en $\vec{v}$. \item (4 pts) Implémenter une fonction Matlab qui résout un système linéaire $A \vec{u} = \vec{f}$ à l'aide de la méthode du gradient. (On peut utiliser le code \texttt{SteepestDescent} dans le polycopi\'e.) R\'esoudre un syst\eme $A \vec{u} = \vec{f}$ o\u $A$ est le Laplacien en dimension deux et utiliser les r\'esultats \a estimer le valeur de $\frac{\kappa(A) - 1}{\kappa(A) + 1}$. \end{enumerate} \end{enumerate} $\;$\\ \\ {\bf \'Evaluation:} \begin{itemize} \item[$\bullet$] Les exercices \item[$\bullet$] Un examen oral durant la session d'examens sur le cours. \end{itemize} La note finale est de : $30 \%$ exercices et $70\%$ examen oral.\\ Les exercices sont notés $15\%$ pour la présentation de séries et $15\%$ pour la présentation au tableau. Il faut passer au tableau au moins 2 fois dans le semestre. \end{document}
 function [u,uk]=SteepestDescent(A,f,u0,tol,maxiter) % STEEPESTDESCENT solves Au=f using the steepest descent method % [u,uk]=SteepestDescent(A,f,u0,tol,maxiter); solves Au=f using steepest % descent starting at the initial guess u0 up to a tolerance tol using % at most maxiter iterations. A has to be symmetric positive definite. % uk contains all the iterates. if nargin<5, maxiter = 100; end if nargin<4, tol =0.0001; end r=f-A*u0; uk(:,1)=u0; k=0; while norm(r)/norm(f)>tol && k
 % Test suite for Serie 3 n = 10; A = rand(n); A = -(A+A')/2; A = diag(-sum(A,2)+0.1) + A; u = 2*rand(n,1)-1; f = A*u; [~,vk]= SteepestDescent(A,f,rand(n,1)); [~,vcol] = size(vk); error = zeros(1,vcol); for i=1:vcol error(i) = norm(vk(:,i) - u); end semilogy(1:vcol,error) %% Laplacian convergence n = 7; A = FDLaplacianTrue(n,2); A = -A; x = linspace(0,1,n); u = cos(x') * sin(x); u = u(:); f = A*u; C = condest(A); C = (C-1)/(C+1); [v,vk] = SteepestDescent(A,f,ones(n^2,1)); [vrow,vcol] = size(vk); error = zeros(vcol,1); for i=1:vcol w = vk(:,i) - u; error(i) = sqrt(w'*A*w); end P = fit((1:vcol)',log(error),'poly1'); plot(1:vcol,log(error),'.-',1:vcol,P(1:vcol),'--') xlabel('Iteration') ylabel('log(Error)') Cest = exp(P.p1); %Cest = (1+Cest)/(1-Cest); \ No newline at end of file
 \documentclass[a4paper]{article} \input{style/head.tex} %------------------------------- % TITLE VARIABLES (identify your work!) %------------------------------- \newcommand{\yourname}{} % replace YOURNAME with your name \newcommand{\yournetid}{} % replace YOURNETID with your NetID \newcommand{\youremail}{} % replace YOUREMAIL with your email \newcommand{\assignmentnumber}{3} % replace X with assignment number \usepackage{preamble} \usepackage{pdfpages} \usepackage{enumitem} \theoremstyle{definition} \newtheorem{ex}{Exercice} \begin{document} %------------------------------- % TITLE SECTION (do not modify unless you really need to) %------------------------------- \input{style/header.tex} %------------------------------- % ASSIGNMENT CONTENT (add your responses) %------------------------------- \begin{ex} \begin{enumerate}[label={\alph*)}] \item Noter que les valeurs propres de Laplacien en dimension deux sont: \begin{equation*} \lambda_{i,j} = -\frac{2}{h^2} \left ( 2 + \cos \left ( \frac{i \pi}{n+1} \right ) + \cos \left ( \frac{j \pi}{n+1} \right ) \right ). \end{equation*} Le plus grand et le plus petit sont $\lambda_{1,1}$ et $\lambda_{n,n}$, respectivement. Rappeler que $\cos(\pi - \theta) = -\cos(\theta)$, alors $\lambda_{1,1} + \lambda_{n,n} = -8/h^2$. Donc, $\alpha_{opt} = -h^2/4$. Consid\'erer le splitting de Laplacien $A=D + L + U$. Dans ce cas, $D^{-1} = \alpha_{opt} I$. Alors, l’itération de Jacobi est \begin{align} u_{k+1} & = -D^{-1} (L+U) u_k + D^{-1} f \\ & = -D^{-1} (A - D) u_k + \alpha_{opt} f \\ & = (I - \alpha_{opt} A) a_k + \alpha_{opt} f \end{align} qui est exactement l'it\'eration de Richardson. \item D'abord, \'ecrire $A=\begin{bmatrix} a & b \\ c & d \end{bmatrix}$. Pour Jacobi: \begin{align} M^{-1} N & = -D^{-1}(U+L) = \begin{bmatrix} 0 & b/a \\ c/d & 0 \end{bmatrix} \\ \implies \rho(M^{-1}N) & = \abs{\sqrt{\frac{bc}{ad}}}. \end{align} On suppose que Jacobi converge, donc $\abs{bc/ad}<1$. Pour GS: \begin{align*} M^{-1} N & = -(D+L)^{-1} U = \begin{bmatrix} 0 & b/a \\ 0 & -bc/ad \end{bmatrix} \\ \implies \rho(M^{-1}N) & = \abs{\frac{bc}{ad}} < 1. \end{align*} \item \end{enumerate} \end{ex} \begin{ex} \begin{enumerate}[label={\alph*)}] \item Consid\'erer la d\'eriv\'ee partielle de $F(u)$ en $u_i$: \begin{align*} \pdxdy{F(u)}{u_i} & = \frac{1}{2} \pdxdy{}{u_i} \begin{bmatrix} u_1 & \dots & u_n \end{bmatrix} A \begin{bmatrix} u_1 \\ \vdots \\ u_n \end{bmatrix} - \pdxdy{}{u_i} \begin{bmatrix} u_1 & \dots & u_n \end{bmatrix} f \\ & = \frac{1}{2} e_i^\top A u + \frac{1}{2} u A e_i^\top - e_i^\top f \\ & = \frac{1}{2} e_i^\top A u + \frac{1}{2} e_i^\top A^\top u - e_i^\top f \\ & = e_i^\top \left( \frac{1}{2} A u + \frac{1}{2} A^\top u - f \right) \\ & = e_i^\top \left( A u - f \right), \end{align*} o\'u $e_i$ est le $i$--i\eme colonne de la matrice idendit\'e. Quand on construit $\nabla F(u)$ avec $\pdxdy{F(u)}{u_i}$ on a $\nabla F(u) = I (Au - f) = Au - f$. \item Soit $\vec{w}$ le vecteur tangent \`a la ligne de niveau de $f(\vec{v})$. La d\'eriv\'ee directionnelle de $f$ dans la direction de $\vec{w}$ est $\nabla f(\vec{v}) \cdot \vec{w}$. Dans cette direction le valeur de $f$ ne change pas, alors $\nabla f(\vec{v}) \cdot \vec{w} = 0$. (Une d\'emonstration irait bien) Soit $\vec{u}$ un vecteur unitaire dans une direction arbitraire qui fait un angle avec $\nabla f(\vec{v})$ de $\theta$. Le produit scalaire entre $\vec{u}$ et $\nabla f(\vec{v})$ est: \begin{align*} \langle \nabla f(\vec{v}), \vec{u} \rangle & = \norm{\nabla f(\vec{v})} \norm{\vec{u}} \cos \theta \\ & = \norm{\nabla f(\vec{v})} \cos \theta \end{align*} qui est minimis\'e par $\theta = -\pi$. Alors, $\vec{u}$ point dans la direction opposite de $\nabla f(\vec{v})$. \end{enumerate} \end{ex} %------------------------------------------------ \end{document}
 % Master Preamble \NeedsTeXFormat{LaTeX2e} \ProvidesPackage{preamble} \RequirePackage{amsmath} \RequirePackage{amssymb} \RequirePackage{amsthm} \RequirePackage{graphicx} \RequirePackage{hyperref} \RequirePackage{subcaption} \newcommand{\dxdy}{\frac{d #1}{d #2}} \newcommand{\dxdyk}{\frac{d^{#3} #1}{d {#2}^{#3}}} \newcommand{\pdxdy}{\frac{\partial #1}{\partial #2}} \newcommand{\liminfty}{\lim_{#1 \to \infty}} \newcommand{\limab}{\lim_{#1 \to #2}} \newcommand{\abs}{\left \vert #1 \right \vert} \newcommand{\norm}{\left \Vert #1 \right \Vert} \newcommand{\order}{\mathcal{O} \left ( #1 \right )} \newcommand{\set}{\left \{ #1 \right \}} \newcommand{\Set}{\left \{ #1 \ \middle \vert \ #2 \right \}} \newcommand{\vmat}{\begin{vmatrix} #1 \end{vmatrix}} \DeclareMathOperator{\sign}{sign} \newcommand{\bbn}{\mathbb{N}} \newcommand{\bbz}{\mathbb{Z}} \newcommand{\bbq}{\mathbb{Q}} \newcommand{\bbr}{\mathbb{R}} \newcommand{\bbc}{\mathbb{C}} \newcommand{\bbf}{\mathbb{F}} \newtheorem{lemma}{Lemma} \newtheorem{thm}{Theorem} \newtheorem{cor}{Corollary} \newtheorem{prop}{Proposition} \newtheorem{conj}{Conjecture} \newtheorem{defn}{Definition} \endinput \ No newline at end of file