// TP 3 //1. Intégration d'une fonction continue sur un intervalle borne function [Imoins, Iplus, N] = Riemann(f, a, b, epsilon) nsubdiv = 20 // Voir la Remarque juste après N = 1 done = %F while (~done) N = 2 * N Imoins = 0 Iplus = 0 for k = 1:N Ik = [ a+(k-1)*(b-a)/N, a+k*(b-a)/N ] // (1) subdiv = linspace(Ik(1), Ik(2), nsubdiv) // Subdivision de l’intervalle fmin = min(f(subdiv)) fmax = max(f(subdiv)) Imoins = Imoins + fmin Iplus = Iplus + fmax end Imoins = Imoins * (b-a)/N Iplus = Iplus * (b-a)/N done = (Iplus - Imoins) < epsilon end // (while) endfunction function y=sinus_pi(x) y=sin(%pi.*x) endfunction [Imoins, Iplus,N] = Riemann(sinus_pi, 0, 1, 10^-4); // Integrale int_0^1 sin(pi t) de = -cos( pi t )/pi ]_0^1= 2/pi z_sin= integrate('sinus_pi(x)', 'x', 0, 1); disp("fonction TP: ", Imoins, Iplus) disp("fonction scilab integrate: ",z_sin) disp("valeur exacte (calculée): ", 2/%pi) // Plot FONCTION sin(pi x) sur 0, 1 // Intégrale de Riemann zx=linspace(0,1,200) zy=sinus_pi(zx) scf() plot(zx,zy) // Méthode des rectangles à gauche avec borne sup // Avec Borne sup (Cf TD) function y=der_sin_pi(x) // fonction dérivée y=%pi.*cos(%pi.*x) endfunction function [I,I_N, dif_I, bor_sup] = Rectangle_g(f, a, b,N) // valeur de l''intégrale, son approximation, leur différence et la borne sup nsubdiv = 20; // deb: borne sup x_sup=linspace(a,b,nsubdiv) y_sup=abs(der_sin_pi(x_sup)) bor_sup=max(y_sup).*((b-a).^2./(2*N)) // fin : borne sup x=linspace(a,b,N) // deb : methode rectangle I=integrate('sinus_pi(x)','x',a,b) I_N = 0; for k = 1:N xk = a+(k-1)*(b-a)/N; I_N = I_N + sinus_pi(xk); end I_N = I_N.*((b-a)./N); dif_I=abs(I_N-I); endfunction [I,I_N, dif_I, bor_sup] = Rectangle_g(sinus_pi, 0, 1,100) // Sans La borne sup (Approche similaire à celle utilisée pour écrire la fonction "Riemann") function [I_rec, I_diff, N] = Rectangle_alt(f, a, b, epsilon) N = 1 done = %F while (~done) N = 2 * N I_rec = 0 for k = 1:N xk= a+(k-1)*(b-a)/N I_rec = I_rec + f(xk) end I_rec= I_rec * (b-a)/N I=integrate('f(x)','x',a,b) I_diff=I - I_rec done = I_diff< epsilon end // (while) endfunction [I_rec, I_diff, N] = Rectangle_alt(sinus_pi, 0, 1, 10^(-4)) // PARTIE 2 CONVERGENCE PONCTUELLE SERIE DE FOURIER FONCTION PERIODIQUE REGULIERE function y=Perio_cont(x) //Graphe de la fonction sur n'importe quel intervalle (plusieurs périodes de longueur 2) y= (abs(x)-(2.*int(abs(x)/2))).*(modulo(floor(abs(x)),2)==0) + (2-(abs(x)-(2.*int(abs(x)/2)))).*(modulo(floor(abs(x)),2)==1) endfunction // Explication //1. fonction paire donc f(-x)=f(x) donc resultat en fonction de abs(x) ; //2. quand x in [ 2k, 2k+1[, k in N alors f(x) = y pour y = x- 2* int(x/2) avec int(x/2)= quotient de la div Eucli de x par 2 ; //condition sur x : ""reste de la div Eu de la partie entiere de x par 2 = 0" et f(x) = x- 2*(quotient de la div euclidienne de x par 2) // 3 quand x in [ 2k+1 , 2(k+1) ] k in N alors f(x) = 2- y (pour se ramener à une valeur entre 0 et 1) // condition sur x : "(modulo(floor(abs(x)),2)==1" scf() zx=linspace(-6,6,601) ; zy=Perio_cont(zx); plot(zx,zy) function y=Perio_cont_alt(x) // Graphe de la fonction uniquement sur [0,2] y= x.*(0<=x & x<1 )+(2-x).*(x>=1 & x<=2) endfunction scf() zx=linspace(0,2,201) ; zy=Perio_cont_alt(zx); plot(zx,zy) // Coefficients de Fourier à la main // a0 =1 // bk = 0 (fonction paire) // ak = -(1-(1)^k)*(2/(k pi)^2) pour k >0 // avec ak = 0 si k>0 pair // ak = -4/((k pi)^2) si k impair // on définit fN = a0/2 + sum_k>=1^N (ak cos( k pi x) + bk sin(k pi x)) // On simplifie alors f_(2N+1)(x)) par 1/2 - 4/(pi^2)sum_{l=0}^N (cos(pi (2l+1) x) )/((2l+1)^2) function y=somme_partielle_Fourier_paire(N,x) dL=0 for l=0:N k=2.*l+1; // le k impair du ak non nul dL=dL+ cos(%pi.*k.*x)/(k.^2); end y=(1/2)-(4/%pi^2).*dL; endfunction zx=linspace(0,2,201) ; zy=Perio_cont(zx); N=[2, 4, 8, 16, 32, 64, 128,256]; zN=zeros(length(N),length(zx)); for i=1:length(N) zN(i,:)=somme_partielle_Fourier_paire(N(i),zx); end zz=[zy', zN']; zx=zx'; scf() plot(zx,zz,LineSpec=1:9) legend(['True','N=2','N=4','N=8','N=16','N=32','N=64','N=128','N=256']) /* Thm (Corollaire de Dirichlet) : les conditions d'application sont vérifiées : f est continue et C1 par morceaux, donc la série de Fourier de f converge normalement vers f sur R et en conséquence elle converge uniformément et simplement vers f sur R */ // 3. "Phénomène de GIBBS POUR UNE FONCTION PERIODIQUE DISCONTINUE " function y=Perio_Discont_supp(x) // uniquement sur le support [-1,1] y= (-1-x).*(-1 <=x & x <0) + (1-x).*(0<=x & x<= 1) endfunction zx=linspace(-1,1,201) ; zy=Perio_Discont_supp(zx) scf() plot(zx,zy) function y=somme_partielle_Fourier_Impaire(N, x) dL=0 for k=1:N dL=dL+2/(k*%pi).*sin(k*%pi.*x) end y=dL endfunction zx=linspace(-1,1,101) ; zy=Perio_Discont_supp(zx); N=[2, 4, 8, 16, 32, 64, 128,256] zN=zeros(length(N),length(zx)); for i=1:length(N) zN(i,:)=somme_partielle_Fourier_Impaire(N(i),zx); end zz=[zy', zN']; zx=zx'; scf() plot(zx,zz,LineSpec=1:9) legend(['True','N=2','N=4','N=8','N=16','N=32','N=64','N=128','N=256'],4) // lower right // Value of fN(0) et ""fN (0) - f(0) " ? zy0N=zeros(length(N),1); dif0N=zeros(length(N),1); for i=1:length(N) zyN0(i)=somme_partielle_Fourier_Impaire(N(i),0); dif0N(i)=zyN0(i)-Perio_Discont_supp(0); end dif0N /* PAS DE CONVERGENCE PONCTUELLE EN x=0 // c''' est cohérent avec le cours (Thm de Dirichlet) en effet f(0)=1 tandis que la série de Fourier de f en x=0 converge vers 0 */ // question 3. Illustration de la cv dans L2 ; on approche la norme L2 au carré de la différence function y=JKN(K,N) y=0 for k=0:K-1 y=y+(abs(Perio_Discont_supp( -1+k*(2/K)) - somme_partielle_Fourier_Impaire(N,-1+k*(2/K))))^2 end y=(2/K)*y; endfunction N=[ 64, 128, 256, 512, 1024, 2048 ] JKN_vect=zeros(length(N),1); K=600 for i=1:length(N) JKN_vect(i,:)=JKN(K,N(i)); end JKN_vect // JKN tend vers 0 avec N qui grandit ce qui est cohérent avec le cours . (cv dans L2)