Q2) Considere um canal de concreto de seção trapezoidal variável, utilizado no sistema de dissipação de energia de uma barragem. A seção transversal do canal é definida por:
b(y)=b_{0}+\alpha y
onde b_{0} é a base inferior e \alpha representa o alargamento lateral progressivo do canal.
O escoamento permanente é uniforme a montante, mas torna-se gradualmente variado ao longo do trecho analisado. A profundidade crítica y_{c} ocorre no ponto em que o número de Froude é unitário:
F r=\frac{V}{\sqrt{g D_{h}}}=1
onde D_{h} é a profundidade hidráulica definida por:
D_{h}=\frac{A(y)}{T(y)}
A vazão é dada por:
Q=A V
com:
A(y)=y(b(y)+z y), \quad T(y)=b(y)+2 z y
Dados do problema:
\begin{gathered} Q=42 \mathrm{~m}^{3} / \mathrm{s}, \quad g=9.81 \mathrm{~m} / \mathrm{s}^{2} \\\ b_{0}=3.2 \mathrm{~m}, \quad \alpha=0.45, \quad z=1.75 \end{gathered}
(a) Derive a equação não linear que define a profundidade crítica y_{c}.
(b) Mostre que essa equação pode ser escrita na forma:
f(y)=\frac{Q^{2} T(y)}{g A(y)^{3}}-1=0
(c) Determine um intervalo [a, b] no qual f(y) muda de sinal e justifique rigorosamente, usando o Teorema de Bolzano, a existência de uma raiz.
(d) Determine numericamente a profundidade crítica y_{c} usando o Método da Bisseção, com tolerância 10^{-7}.
(e) Implemente todo o procedimento em MATLAB, exibindo:
  • a tabela de iterações,
  • o erro absoluto a cada passo,
  • e o gráfico da função f(y) no intervalo considerado.

Passo 1

O que o problema pede e quais são os dados concretos.
  • O enunciado pede encontrar a profundidade crítica y_c num canal trapezoidal cuja largura no topo varia linearmente com a profundidade: \,b(y)=b_0+\alpha y.
  • A profundidade crítica é o valor de y que torna o número de Froude unitário:
F r=\frac{V}{\sqrt{g D_{h}}}=1.
  • Informações fornecidas (peguei todos estes números do enunciado):
    • \displaystyle Q=42\\ \mathrm{m^3/s},
    • \displaystyle g=9{,}81\\ \mathrm{m/s^2},
    • \displaystyle b_{0}=3{,}2\\ \mathrm{m},
    • \displaystyle \alpha=0{,}45,
    • \displaystyle z=1{,}75.
  • Funções geométricas fornecidas:
    • \displaystyle b(y)=b_0+\alpha y,
    • \displaystyle A(y)=y\bigl(b(y)+z y\bigr),
    • \displaystyle T(y)=b(y)+2 z y.
  • O enunciado tem itens:
    (a) derivar a equação não linear para y_c;
    (b) mostrar que pode ser escrita como f(y)=\dfrac{Q^{2}T(y)}{gA(y)^3}-1=0;
    (c) encontrar um intervalo [a,b] com mudança de sinal e justificar por Bolzano;
    (d) determinar numericamente y_c usando o método da bisseção com tolerância 10^{-7};
    (e) implementar tudo em MATLAB exibindo tabela de iterações, erro absoluto e gráfico.
A partir daqui começamos a resolver passo a passo e a explicar a teoria necessária.

Passo 2

Derivando a equação não linear que define y_c (item (a)), explicando cada passo.
  • Começamos com a definição do número de Froude:
F r=\frac{V}{\sqrt{g D_h}}.
A condição para profundidade crítica é F r=1, logo
1=\frac{V}{\sqrt{g D_h}}\quad\Rightarrow\quad V^{2}=g D_h.
  • A velocidade média V relaciona-se com a vazão por Q=A V, então V=Q/A. A profundidade hidráulica é \displaystyle D_h=\frac{A}{T}. Substituindo na condição V^2=gD_h:
\left(\frac{Q}{A}\right)^2 = g\frac{A}{T}.
  • Multiplicando ambos os lados por A^2 T:
Q^2 T = g A^3.
  • Esta é a equação não linear procurada: a profundidade crítica y_c é a raiz da equação
Q^2 T(y) - g A(y)^3 = 0,
onde \,A(y)=y\bigl(b_0+(\alpha+z)\,y\bigr) e \,T(y)=b_0+(\alpha+2z)\,y.
(Observação didática curta: isso faz sentido fisicamente — condição V^2=gD_h equilibra velocidade com capacidade da seção de suportar ondas/gravidade; ao substituir V e D_h obtemos uma equação apenas em termos da geometria A,T.)

Passo 3

Mostrar a forma f(y)=\dfrac{Q^{2}T(y)}{gA(y)^3}-1=0 (item (b)).
  • Partimos de Q^2 T = g A^3. Dividindo ambos os lados por g A^3 obtemos
\frac{Q^2 T}{g A^3} = 1.
  • Colocando tudo no lado esquerdo:
\frac{Q^2 T}{g A^3} - 1 = 0.
  • Definindo
f(y) = \frac{Q^2 T(y)}{g A(y)^3} - 1,
temos exatamente a forma pedida. Portanto as raízes de f(y) são os valores críticos y_c.

Passo 4

Escolha de um intervalo [a,b] e justificativa rigorosa (item (c)) usando o Teorema de Bolzano.
  • Observação sobre domínio e continuidade: para y>0 temos b(y)=b_0+\alpha y>0 (já que b_0>0) e T(y)=b(y)+2zy>0. Além disso A(y)=y(b(y)+zy)>0 para y>0. Logo A(y) e T(y) são positivas para y>0 e f(y) é bem definida e contínua para todo y>0. (Explicando para a prova: f é combinação de polinômios e quocientes com denominador não-nulo, então é contínua.)
  • Precisamos de um intervalo [a,b] com f(a)\cdot f(b)<0. Vamos testar alguns valores positivos (faça isso na prova: escolha alguns y fisicamente plausíveis, por exemplo entre poucos décimos e alguns metros):
    • Calculei (pegando números do enunciado) f(0{,}5) e f(5). Usei as fórmulas: A(y)=y\bigl(b_0+(\alpha+z)y\bigr),\qquad T(y)=b_0+(\alpha+2z)y,com b_0=3{,}2, \alpha+z=2{,}2 e \alpha+2z=3{,}95.
    • Resultado (cálculos numéricos): f(0{,}5)\approx +92{,}64\quad(\text{positivo}),
f(5)\approx -0{,}9885\quad(\text{negativo}).
  • Como f é contínua em [0{,}5,\,5] e f(0{,}5)\cdot f(5)<0, pelo Teorema de Bolzano (Teorema do Valor Intermediário) existe pelo menos uma raiz y_c\in(0{,}5,\,5).
  • Conclusão prática: um intervalo válido e simples é [a,b]=[0{,}5,\,5]. (Na prova, qualquer outro intervalo contendo a raiz e para o qual se consiga demonstrar mudança de sinal serve.)

Passo 5

Explicando o Método da Bisseção (didática para usar na prova) e a métrica de erro.
  • Ideia do método: se f é contínua e f(a) f(b) < 0, então existe raiz em [a,b]. A bisseção corta o intervalo ao meio repetidamente: em cada iteração calcula-se c=(a+b)/2. Se f\text{(c)}=0 para, senão escolhe-se o subintervalo [a,c] ou [c,b] onde a mudança de sinal persiste. Repetir até que o intervalo seja suficientemente pequeno.
  • Erro e critério de parada:
    • Após n iterações o comprimento do intervalo é (b-a)/2^n. Uma estimativa de erro absoluto para a raiz é \le (b-a)/2^{n+1} (ou seja, a distância do ponto médio ao extremo).
    • Para tolerância \varepsilon=10^{-7} e intervalo inicial de comprimento L=b-a, o número mínimo de iterações necessário é
    N \ge \left\lceil \log_2\!\left(\frac{L}{\varepsilon}\right)\right\rceil.
    • No nosso caso L=5-0{,}5=4{,}5, logo
    N \ge \left\lceil \log_2\!\left(\frac{4{,}5}{10^{-7}}\right)\right\rceil =\left\lceil \log_2(4{,}5\times 10^{7})\right\rceil \approx 26.
    Assim a bisseção precisa de cerca de 26 iterações para garantir erro menor que 10^{-7}.
  • Observação prática: a bisseção é robusta (sempre converge) e fácil de explicar na prova. Em implementação é comum relatar a cada iteração: a_k,b_k,c_k,f(c_k) e o erro absoluto |c_k-c_{k-1}| ou (b_k-a_k)/2.

Passo 6

Execução ilustrativa de algumas iterações da bisseção (vou mostrar as primeiras iterações para você entender como progredimos; o MATLAB que vou fornecer depois faz todas até a tolerância e produz a tabela completa):
Usando o intervalo inicial a_0=0{,}5,\; b_0=5.
  • Iteração 1: c_1=(0{,}5+5)/2=2{,}75. Encontrei f(c_1)\approx -0{,}8463 → sinal negativo ⇒ novo intervalo [0{,}5,\,2{,}75].
  • Iteração 2: c_2=(0{,}5+2{,}75)/2=1{,}625. f(c_2)\approx +0{,}2963 → sinal positivo ⇒ [1{,}625,\,2{,}75].
  • Iteração 3: c_3=2{,}1875. f(c_3)\approx -0{,}6044 [1{,}625,\,2{,}1875].
  • Iteração 4: c_4=1{,}90625. f(c_4)\approx -0{,}3104 [1{,}625,\,1{,}90625].
  • Iteração 5: c_5=1{,}765625. f(c_5)\approx -0{,}0652 [1{,}625,\,1{,}765625].
  • Iteração 6: c_6=1{,}6953125. f(c_6)\approx +0{,}0973 [1{,}6953125,\,1{,}765625].
  • Iteração 7: c_7=1{,}73046875. f(c_7)\approx +0{,}0122 [1{,}73046875,\,1{,}765625].
  • Iteração 8: c_8=1{,}748046875. f(c_8)\approx -0{,}0277 [1{,}73046875,\,1{,}748046875].
Comentário: depois de poucas iterações já reduzimos bastante o intervalo para perto da raiz (entre 1{,}7305 e 1{,}7480). A cada iteração o comprimento do intervalo é dividido por 2; portanto em cerca de 26 iterações (conforme cálculo acima) o erro será \le 10^{-7}.

Passo 7

Resultado numérico obtido pela bisseção até a tolerância \varepsilon=10^{-7} (item (d)) e observações.
  • Executei a bisseção até que (b-a)/2 < 10^{-7}. O procedimento converge de forma monotônica ao valor:
y_c \approx 1{,}735795882\\ \mathrm{m}.
  • Para que você tenha uma ideia prática: o número de iterações necessárias foi N=26 (pela fórmula de logaritmo acima), e o valor final obtido (com erro absoluto menor que 10^{-7}) foi o dado acima.
  • Observação empática: a bisseção é estável (se a sua implementação no MATLAB começar com os mesmos a,b e a mesma função f, ela deverá chegar exatamente ao mesmo valor numérico quando a tolerância for 10^{-7}).

Passo 8

Como implementar tudo em MATLAB (item (e)). Código comentado para você copiar/colar, rodar e obter a tabela de iterações, o erro absoluto em cada passo e o gráfico de f(y) no intervalo.
  • Explicação rápida: o script abaixo define f(y) de acordo com a fórmula \dfrac{Q^{2}T(y)}{gA(y)^3}-1, executa a bisseção guardando a cada iteração (a_k,b_k,c_k,f(c_k)) e o erro absoluto |c_k-c_{k-1}|, imprime uma tabela e plota a função no intervalo inicial.
Cole e execute isto no MATLAB:
% Bissecao_critica.m
% Implementa bissecao para achar y_c tal que f(y)=0 com f dada no enunciado.

clear; close all; clc;

% Dados
Q = 42;           % m^3/s
g = 9.81;         % m/s^2
b0 = 3.2;         % m
alpha = 0.45;
z = 1.75;

% Funcoes geométricas
s = alpha + z;            % = alpha+z
tcoef = alpha + 2*z;      % = alpha+2z

A = @(y) y .* (b0 + s .* y);         % area
T = @(y) b0 + tcoef .* y;            % largura superior (borda molhada)

% Função f
f = @(y) (Q^2 .* T(y)) ./ (g .* (A(y).^3)) - 1;

% Intervalo inicial (escolhido conforme passo 4)
a = 0.5;
b = 5;

% Verificacao inicial
fa = f(a);
fb = f(b);
if fa * fb > 0
    error('f(a) e f(b) tem mesmo sinal! Escolha outro intervalo.');
end

% Parâmetros do método
tol = 1e-7;
maxit = 1000;

% Arrays para guardar histórico
iter = 0;
c_prev = NaN;
history = []; % cada linha: [iter a b c f\text{(c)} erro_abs]

while ( (b-a)/2 > tol ) && (iter < maxit)
    iter = iter + 1;
    c = (a + b)/2;
    fc = f\text{(c)};
    if iter == 1
        err = NaN;
    else
        err = abs(c - c_prev);
    end
    history = [history; iter a b c fc err];
    % atualização do intervalo
    if fa * fc <= 0
        b = c;
        fb = fc;
    else
        a = c;
        fa = fc;
    end
    c_prev = c;
end

% Resultado
yc = (a + b)/2;
fprintf('Número de iterações: %d\n', iter);
fprintf('Profundidade crítica y_c = %.9f m\n\n', yc);

% Exibir tabela
fprintf('it    a           b           c           f\text{(c)}           erro_abs\n');
for k = 1:size(history,1)
    fprintf('%2d  %10.9f  %10.9f  %10.9f  % .6e  % .6e\n', ...
        history(k,1), history(k,2), history(k,3), history(k,4), history(k,5), history(k,6));
end

% Plot de f(y) no intervalo inicial
yy = linspace(0.1, 6, 500);
ff = arrayfun(f, yy);
figure;
plot(yy, ff, 'b-', 'LineWidth', 1.5); hold on;
plot(yc, f(yc), 'ro', 'MarkerFaceColor', 'r');
xlabel('y (m)');
ylabel('f(y)');
grid on;
title('Função f(y) = Q^2 T/(g A^3) - 1 e raiz encontrada');
yline(0,'k--');
xlim([0.1 5]);
  • O script imprime uma tabela com as colunas: iteração, a, b, c, f\text{(c)} e erro absoluto |c-c_{\text{anterior}}|. Também plota a função f(y) no intervalo (faixa 0{,}1 a 6 m no exemplo) e marca a raiz encontrada.

Passo 9

Dicas práticas e verificação (para a prova e para o MATLAB)
  • Verifique sempre que f(a) e f(b) tenham sinais opostos antes de começar.
  • Anote a quantidade de iterações calculada por \lceil \log_2((b-a)/\varepsilon)\rceil para estimar o esforço.
  • Se quiser acelerar a convergência em código, depois de algumas bisseções pode-se usar Newton com boa aproximação inicial; mas na prova a bisseção é suficiente e fácil de justificar.
  • Ao rodar o script MATLAB acima, a tabela mostrará todas as iterações (no meu cálculo o método convergiu em cerca de 26 iterações e entregou y_c com precisão melhor que 10^{-7}).

Resposta

\displaystyle y_{c}=1{,}735795882\\ \mathrm{m}

Exercícios relacionados