Análise Hidrodinâmica do Navio I

...

psi(i,j) = -3.6*y(i,j) -2/pi*(atan(y(i,j)/(x(i,j) + 4))-pi) -2/pi*(atan(y(i,j)/(x(i,j) - 3))-pi) +5/pi*(atan(y(i,j)/(x(i,j) - 5))-pi);

end

end

end

end

end

else

psi(i,j) = -3.6*y(i,j) -2/pi*atan(y(i,j)/(x(i,j) + 4)) -2/pi*atan(y(i,j)/(x(i,j) - 3)) +5/pi*atan(y(i,j)/(x(i,j) - 5));

end

end

end

Esta sequência de if’s garante que os ângulos sejam computados corretamente no código, pois como a função arco-tangente retorna sempre valores no intervalo , foi necessário ajustar os resultados da função para representar o real valor de , medido no sentido anti-horário em relação ao semi-eixo positivo de . Dessa forma, o plano cartesiano foi dividido em diferentes regiões, definidas pela localização da fonte e dos sumidouros e pelos semi-planos positivo e negativo. Dependendo da região de análise, a função arco-tangente poderia ser ajustada em ou , resultando no valor de .[pic 13][pic 14][pic 15][pic 16][pic 17][pic 18][pic 19][pic 20][pic 21]

[vy,vx]=gradient(psi,.03);

vy=-vy;

V=sqrt(vx.^2+vy.^2);

P= 1/2*1000*(3.6^2-(vx.^2+vy.^2));

Essas linhas correspondem à utilização do gradiente da função corrente para descobrir e e à implementação da equação de Bernoulli aplicada ao escoamento, como apresentada anteriormente.[pic 22][pic 23][pic 24]

m=zeros(fix((size(x)-1)/40)+1);

n=zeros(fix((size(x)-1)/40)+1);

u=zeros(fix((size(x)-1)/40)+1);

v=zeros(fix((size(x)-1)/40)+1);

for i=0:40:size(x)

for j=0:40:size(x)

m(i/40+1,j/40+1)=x(i+1,j+1);

n(i/40+1,j/40+1)=y(i+1,j+1);

u(i/40+1,j/40+1)=vx(i+1,j+1);

v(i/40+1,j/40+1)=vy(i+1,j+1);

end

end

Esta etapa foi necessária para reduzir a quantidade de vetores a serem representados mais adiante no código na etapa de plotagem do campo de velocidades, ou seria plotado um vetor para cada ponto em que a velocidade foi calculada pelo código. Para evitar isso, foram criados vetores contendo informações de apenas alguns pontos e suas respectivas velocidades.

figure('name','Linhas de corrente','numbertitle', 'off')

contour(x,y,psi,40)

hold on

quiver(m,n,u,v,40)

figure('name','Distribuição de pressões','numbertitle', 'off')

c=-7000:20:6500;

contour(x,y,P,c)

hold on

h=[0 0];

contour(x,y,psi,h)

Resultados, comentários

Abaixo, são apresentadas as imagens que representam graficamente os resultados encontrados.

[pic 25]

[pic 26]

Na primeira, podemos observar as linhas de corrente do escoamento com o campo de velocidades. Na segunda, a distribuição de pressões no escoamento, com a representação da linha de corrente , que pode ser utilizada para delimitar o formato de um corpo sólido no escoamento.[pic 27]

Para informações mais detalhadas, também é possível consultar os valores numéricos das matrizes resultantes ao final da execução do código.

É possível observar pontos de descontinuidade na fonte e nos sumidouros, perto dos quais as velocidades tendem a infinito (e as pressões a menos infinito). Fisicamente, isso pode ser explicado pela modelagem da fonte e dos sumidouros como pontos. Por essa hipótese, certa quantidade de fluido converge para um único ponto geométrico, o que exige que essa quantidade de fluido atinja velocidades cada vez mais altas quanto mais ela se aproximasse desse ponto. Considerando a análise na qual a linha de corrente zero é o contorno de um corpo sólido, esse problema seria suprimido, pois a fonte e os sumidouros ficariam na região interior a essa linha de