Trabalho 02 - Cálculo Numérico

Prévia do material em texto

<p>UNIVERSIDADE FEDERAL DO RIO DE JANEIRO</p><p>CAMPUS MACAÉ PROF. ALOISIO TEIXEIRA</p><p>CÁLCULO NUMÉRICO</p><p>PROF. FELIPE SALES</p><p>TRABALHO 02: CÁLCULO NUMÉRICO</p><p>ANDRÉ FELIPE DE OLIVEIRA CARVALHO</p><p>DRE: 118024140</p><p>MACAÉ - RIO DE JANEIRO</p><p>2019</p><p>Projeto 4.4</p><p>O problema de se determinar um polinômio:</p><p>,</p><p>de grau no máximo , tal que:</p><p>,</p><p>onde são constantes, pode ser resolvido através da obtenção da solução de um sistema linear. Determine um polinômio de grau 3, que satisfaça a condição acima, considerando , e:</p><p>A.</p><p>B.</p><p>Resolvendo o sistema linear resultante por método numérico à sua escolha.</p><p>Como o que se quer é um polinômio de no máximo terceiro grau, ou seja, , temos que .</p><p>Letra A:</p><p>R = .</p><p>Para temos que , ou seja,</p><p>Resolvendo a integral temos: .</p><p>Analogamente,</p><p>Para temos ;</p><p>Para temos ;</p><p>Para temos .</p><p>Resolvendo o sistema utilizando-se do método de Newton, obtemos , , e .</p><p>Logo, temos que o polinômio é .</p><p>Letra B:</p><p>R = .</p><p>Para temos que , ou seja, .</p><p>Resolvendo a integral temos: .</p><p>Analogamente,</p><p>Para temos ;</p><p>Para temos ;</p><p>Para temos .</p><p>Resolvendo o sistema utilizando-se do método de Newton, obtemos , , e .</p><p>Logo, temos que o polinômio é aproximadamente .</p><p>Algoritmo implementado no Scilab:</p><p>Para a Letra A:</p><p>y(1)=2*x(1)+(2/3)*x(3) - 2/3</p><p>y(2)=(2/3)*x(2)+(2/5)*x(4) - 4/3</p><p>y(3)=(2/3)*x(1)+(2/5)*x(3) - 6/5</p><p>y(4)=(2/5)*x(2)+(2/7)*x(4) - 4/5</p><p>Para a Letra B:</p><p>y(1)=2*x(1)+(2/3)*x(3) - 2</p><p>y(2)=(2/3)*x(2)+(2/5)*x(4) - 2</p><p>y(3)=(2/3)*x(1)+(2/5)*x(3) - 2/3</p><p>y(4)=(2/5)*x(2)+(2/7)*x(4) - 58/35</p><p>//Método de Newton para sistemas não lineares</p><p>clear; clc</p><p>//Entrada</p><p>function y=F(x)</p><p>y(1)=...</p><p>y(2)=...</p><p>y(3)=...</p><p>y(4)=...</p><p>endfunction</p><p>function y=J(x) //Matriz Jacobiana</p><p>y = [2 0 2/3 0;0 2/3 0 2/5;</p><p>2/3 0 2/5 0;0 2/5 0 2/7]</p><p>endfunction</p><p>erro = 0.001</p><p>x0 = [1 1 1 1]' //Palpite</p><p>pare = %f</p><p>x = x0</p><p>n = length(F(x0))</p><p>while pare == %f do</p><p>//Resolver o sistema linear</p><p>A = [J(x) -F(x)]</p><p>for i =1:n do</p><p>m = i //Pivotamento parcial</p><p>for k =(i+1):n do</p><p>if abs(A(k,i)) > abs(A(m,i)) then</p><p>m = k</p><p>end</p><p>end</p><p>L = A(i,:) //Troca de linhas</p><p>A(i,:) = A(m,:)</p><p>A(m,:) = L</p><p>A(i,:) = A(i,:)/A(i,i)</p><p>for k =1:n do</p><p>if k ~=i then</p><p>A(k,:) = A(k,:) - (A(k,i)/A(i,i))*A(i,:)</p><p>end</p><p>end</p><p>end</p><p>r = A(:,n+1)</p><p>y = x + r</p><p>if max(abs(y-x))<erro*max(abs(y)) then</p><p>pare = %t</p><p>end</p><p>x=y</p><p>end</p><p>disp("A solução é:")</p><p>disp(y)</p><p>Projeto 5.3</p><p>Suponha uma barra de metal homogêneo, como na figura a seguir, onde AB = CD = 4 : AC = BD = 3.</p><p>A temperatura ao longo de AB, AC, BD é mantida constante e igual a 0°C, enquanto que ao longo de CD ela é igual a 1°C. A distribuição do calor na barra obedece à seguinte equação:</p><p>com as condições de contorno:</p><p>para ,</p><p>para ,</p><p>para ,</p><p>para .</p><p>A solução numérica deste problema pode ser obtida considerando-se uma divisão do retângulo ABCD em retângulos menores a partir de uma divisão de AB em intervalos iguais de amplitude e de uma divisão de CD em intervalos iguais de amplitude , como é mostrado na figura a seguir:</p><p>Nessa figura estamos considerando . A temperatura nos pontos internos pode ser obtida numericamente simulando as derivadas segundas de (5.25), pelas diferenças de segunda ordem ∆2 de</p><p>modo que para , obtemos:</p><p>+ = 0</p><p>para cada par em . Assim, por exemplo, para o ponto da figura anterior vale:</p><p>+ = 0 .</p><p>Considerando todos os pontos da figura anterior obtemos um sistema de 6 equações lineares nas incógnitas: . Utilizando a fórmula da temperatura nos pontos internos, obtém-se as seguintes equações:</p><p>;	 ;</p><p>;		 ;</p><p>;	.</p><p>Utilizando o método de Newton para resolver o sistema, obtém-se as temperaturas (em ºC) nos seguintes pontos:</p><p>A solução é:</p><p>0.7971014</p><p>0.6521739</p><p>0.4637681</p><p>0.5362319</p><p>0.3478261</p><p>0.2028986</p><p>--> F(y)</p><p>ans =</p><p>-4.441D-16</p><p>0.</p><p>0.</p><p>-2.220D-16</p><p>0.</p><p>0.</p><p>Algoritmo implementado no Scilab</p><p>//Método de Newton para sistemas não lineares</p><p>clear; clc</p><p>function y=F(x)</p><p>y(1)=-4*x(1)+x(2)+x(4)+2</p><p>y(2)=x(1)-4*x(2)+x(3)+x(5)+1</p><p>y(3)=x(2)-4*x(3)+x(6)+1</p><p>y(4)=x(1)-4*x(4)+x(5)+1</p><p>y(5)=x(2)+x(4)-4*x(5)+x(6)</p><p>y(6)=x(3)+x(5)-4*x(6)</p><p>endfunction</p><p>function y=J(x) //Matriz Jacobiana</p><p>y=[-4 1 0 1 0 0;1 -4 1 0 1 0;0 1 -4 0 0 1;1 0 0 -4 1 0;0 1 0 1 -4 1;0 0 1 0 1 -4]</p><p>endfunction</p><p>erro = 0.000001</p><p>x0 = [1 1 1 1 1 1]' //Palpite</p><p>pare = %f</p><p>x = x0</p><p>n = length(F(x0))</p><p>while pare == %f do</p><p>A = [J(x) -F(x)]</p><p>for i =1:n do</p><p>m = i //Pivotamento parcial</p><p>for k =(i+1):n do</p><p>if abs(A(k,i)) > abs(A(m,i)) then</p><p>m = k</p><p>end</p><p>end</p><p>L = A(i,:) //Troca de linhas</p><p>A(i,:) = A(m,:)</p><p>A(m,:) = L</p><p>A(i,:) = A(i,:)/A(i,i)</p><p>for k =1:n do</p><p>if k ~=i then</p><p>A(k,:) = A(k,:) - (A(k,i)/A(i,i))*A(i,:)</p><p>end</p><p>end</p><p>end</p><p>r = A(:,n+1)</p><p>y = x + r</p><p>if max(abs(y-x))<erro*max(abs(y)) then</p><p>pare = %t</p><p>end</p><p>x=y</p><p>end</p><p>disp("A solução é:")</p><p>disp(y)</p><p>Anexo 2A</p><p>Considere o sistema abaixo. Determine se existe alguma solução do sistema nas limitações indicadas. Se existir, calcule pelo menos uma solução.</p><p>; ; .</p><p>A solução é:</p><p>0.3333765</p><p>-0.0152686</p><p>-0.5238506</p><p>--> F(y)</p><p>ans =</p><p>0.0001616</p><p>-0.7316834</p><p>0.000066</p><p>Algoritmo implementado no Scilab:</p><p>//Método de Newton para sistemas não lineares</p><p>clear;clc</p><p>function y=F(x)</p><p>y(1)= 3*x(1)-cos(x(2)*x(3))-0,5</p><p>y(2)= 4*x(1)**2-625*x(2)**2+2*x(2)-1</p><p>y(3)= exp(-x(1)*x(2))+20*x(3)+((10*%pi-3)/3)</p><p>endfunction</p><p>function y=J(x) //Matriz Jacobiana</p><p>y = [3,sin(x(2)*x(3))*x(3),sin(x(2)*x(3))*x(2);</p><p>8*x(1),-(1250*x(2))+2,0;</p><p>exp(-x(1)*x(2))*(-x(2)),exp(-x(1)*x(2))*(-x(1)),20]</p><p>endfunction</p><p>erro = 0.000001</p><p>x0 = [1,1,1]'</p><p>pare = %f</p><p>x = x0</p><p>c = 0</p><p>n = length(F(x0))</p><p>while pare == %f do</p><p>A = [J(x) -F)(x]</p><p>for i =1:n do</p><p>m = i</p><p>for k =(i+1):n do</p><p>if abs(A(k,i)) > abs(A(m,i)) then</p><p>m = k</p><p>end</p><p>end</p><p>L = A(i,:)</p><p>A(i,:) = A(m,:)</p><p>A(m,:) = L</p><p>A(i,:) = A(i,:)/A(i,i)</p><p>for k =1:n do</p><p>if k ~=i then</p><p>A(k,:) = A(k,:) - (A(k,i)/A(i,i))*A(i,:)</p><p>end</p><p>end</p><p>end</p><p>r = A(:,n+1)</p><p>y = x + r</p><p>if max(abs(y-x))<erro*max(abs(y)) then</p><p>pare = %t</p><p>end</p><p>if c >100 then</p><p>pare = %t</p><p>end</p><p>x=y</p><p>c = c+1</p><p>end</p><p>disp("A solução é:")</p><p>disp(y)</p><p>Anexo 1D</p><p>Considere o sistema abaixo. Determine se existe alguma solução do sistema. Se existir, calcule-a. Se existir mais do que uma, calcule pelo menos duas soluções.</p><p>Com x0 = [10 10]':</p><p>A solução é:</p><p>0.121242</p><p>0.2711052</p><p>--> F(y)</p><p>ans =</p><p>0.0000001</p><p>5.018D-10</p><p>Com x0 = [-100 -100]':</p><p>A solução é:</p><p>-0.0981634</p><p>0.2195001</p><p>--> F(y)</p><p>ans =</p><p>1.469D-11</p><p>-1.601D-14</p><p>Algoritmo implementado no Scilab:</p><p>//Método de Newton para sistemas não lineares</p><p>clear; clc</p><p>function y=F(x)</p><p>y(1)= 5*x(1)**2-x(2)**2</p><p>y(2)= x(2)-0.25*(sin(x(1))+cos(x(2)))</p><p>endfunction</p><p>function y=J(x)</p><p>y = [10*x(1) -2*x(2);</p><p>-0.25*cos(x(1)) 1+0.25*sin(x(2))]</p><p>endfunction</p><p>erro = 0.001</p><p>x0 = ...</p><p>pare = %f</p><p>x = x0</p><p>c = 0</p><p>n = length(F(x0))</p><p>while pare == %f do</p><p>A = [J(x) -F(x)]</p><p>for i =1:n do</p><p>m = i</p><p>for k =(i+1):n do</p><p>if abs(A(k,i)) > abs(A(m,i)) then</p><p>m = k</p><p>end</p><p>end</p><p>L = A(i,:)</p><p>A(i,:) = A(m,:)</p><p>A(m,:) = L</p><p>A(i,:) = A(i,:)/A(i,i)</p><p>for k =1:n do</p><p>if k ~=i then</p><p>A(k,:) = A(k,:) - (A(k,i)/A(i,i))*A(i,:)</p><p>end</p><p>end</p><p>end</p><p>r = A(:,n+1)</p><p>y = x + r</p><p>if max(abs(y-x))<erro*max(abs(y)) then</p><p>pare = %t</p><p>end</p><p>if c >100 then</p><p>pare = %t</p><p>end</p><p>x=y</p><p>c = c+1</p><p>end</p><p>disp("A solução é:")</p><p>disp(y)</p><p>image3.png</p><p>image5.png</p><p>image4.png</p><p>image6.png</p>