Diferencia entre revisiones de «Visualización y análisis del comportamiento físico de una placa sometida a distintos campos (Grupo 16-C)»

Revisión del 21:14 6 dic 2013

Este artículo está en versión beta. El autor de este artículo no lo ha terminado todavía, por favor no lo edites hasta que elimine este mensaje.

Mallado de una placa en forma de anillo

1 Introducción

Nuestro estudio se ha centrado en el análisis del comportamiento de un sólido plano con forma de corona circular. Para ello definiremos los distintos campos que afectan a dicho cuerpo, y para facilitar la comprensión de su comportamiento físico, representaremos gráficamente estos campos. Con el fin de realizar estas representaciones nos serviremos del programa informático MatLAB.

2 Superficie de estudio

La superficie sobre la que hemos realizado el estudio es una placa plana con forma de corona circular, cuyos radios exterior e interior son respectivamente 2 y 1. La representación de esta placa se ha hecho en MatLAB mediante el siguiente código:

% definición de parámetros de la superficie
u = [1:0.1:2];          
v = [0:0.1:2*pi+0.1];
%creación del mallado
[uu,vv]=meshgrid(u,v); 
% parametrización
xx=uu.*cos(vv);        
yy=uu.*sin(vv);
mesh(xx,yy,0*xx)      
axis([-3,3,-3,3]); view(2);

3 Campo de temperaturas

Partimos del campo de temperaturas [math]T(x,y)=e^{-y}[/math] . Un campo escalar de temperaturas depende de la posición de cada punto, es decir de sus coordenadas.

3.1 Distribución de temperaturas

El campo de temperaturas depende en nuestro caso en exclusiva de la variable y, es decir, el valor de la temperatura es independiente al valor de la coordenada x. De este hecho, podemos deducir que una vez fijemos un valor para la coordenada y en la superficie de estudio, la temperatura se mantendrá constante para cualquier valor de x. Sabiendo esto podríamos hacernos una idea de cómo sería la distribución de las temperaturas a lo largo de la placa. Además podemos observar fácilmente que al tratarse de una función exponencial, la separación de las curvas de nivel irá variando de una manera geométrica. La representación del campo de temperaturas en MatLAB nos permitirá confirmar estas afirmaciones y nos dará una visión más clara de la distribución de temperaturas.

T = e.^(-yy); 	%% Campo escalar temperatura
subplot(1,2,1);
mesh(xx,yy,T);    % Parametrización
view(2)
axis([-3,3,-3,3])  
subplot(1,2,2);
mesh(xx,yy,T); colorbar;  %visualización 3D

3.2 Variación de temperaturas

Para estudiar la variación de temperaturas a lo largo de la placa hay que analizar el gradiente [math]\nabla T[/math], que en función del campo de temperaturas dado, será: [math]\nabla T = -e^{-y} [/math] Una vez calculado este campo vectorial, podemos proceder a representarlo sobre la superficie, conjuntamente con las curvas de nivel:

Gradiente= -e.^(-yy);	%% Gradiente del campo

subplot(1,2,1)
hold on
contour(xx,yy,T,20)
mesh(xx,yy,T);
quiver3(xx,yy,T,0*xx,Gradiente,0*T);
hold off

subplot(1,2,2)
hold on
quiver(xx,yy,0*xx,Gradiente); colorbar
contour(xx,yy,T,20)
hold off

De la visualización extraemos dos conclusiones. En primer lugar, observamos que el módulo del gradiente va aumentando conforme nos desplazamos en sentido negativo a lo largo del eje de ordenadas. En segundo lugar y más importante, se observa la perpendicularidad entre las curvas de nivel y el gradiente del campo de temperaturas. La razón de esta perpendicularidad radica en la propia definición de ambos conceptos: las curvas de nivel son el lugar geométrico de puntos equipotenciales, es decir, que el campo les asigna igual valor de temperatura, y el gradiente indica la dirección y sentido del incremento de los valores de dicho campo en cada punto. Por lo tanto, el gradiente no puede tener componente paralelo a la curva de nivel, pues esto supondría que el campo crece en tal dirección, lo que contradice las definiciones previas.

4 Campo de desplazamiento

Sobre la superficie de estudio aplicamos una percusión que produce sobre ella un desplazamiento en un instante [math]t_0[/math] que quedará definido por el vector:

[math]
\vec u(\rho,\theta)=\frac{\sin(\pi \theta)}{20\rho^2}\vec g_{\theta}.
[/math]

Aplicando este campo a cada punto obtendremos su desplazamiento respecto a la posición original.

4.1 Consecuencias de la percusión

El resultado de la percusión sobre el sólido es un desplazamiento definido por el vector [math]\vec u(\rho,\theta) [/math]. Representando en MatLAB este campo vectorial sobre la placa, podremos observar en qué dirección se desplazará cada punto. Implementando el siguiente código:

i=1; j = 1;	
Ux = [];
Uy = [];

for i= 1:64	% bucle para calcular el vector de cada punto
  for j = 1:11
    Ux(i,j)=(-sin(pi*vv(i,j)).*sin(vv(i,j)))./(20.*uu(i,j));
    Uy(i,j)=((sin(pi.*vv(i,j))).*cos(vv(i,j)))./(20.*uu(i,j));
    j = j+1;
    end
  i= i+1;
end
quiver(xx,yy,Ux,Uy);	% Representación del campo vectorial

A primera vista, podemos apreciar que el módulo de cada vector [math]\vec u[/math] en cada punto, es muy pequeño, por lo que es complicado deducir cuál será el comportamiento del sólido debido a la percusión. Con el fin de facilitar la comprensión de los vectores, multiplicamos al vector [math]\vec u[/math] por diez, así estudiaremos una percusión de mayor intensidad en la que los desplazamientos en el sólido serán más notables.

Con este nuevo campo, apreciamos que se produce un desplazamiento de materia hacia tres zonas radiales, así como una deformación del contorno de la placa. Recuperando el módulo original de [math]\vec u[/math], aplicamos dicho vector al sólido para obtener la posición de la placa tras la percusión.

subplot(1,2,1); 		%% situación inicial
mesh(xx,yy,0*xx);
axis([-3,3,-3,3])
view(2)

i=1; j = 1; 		%% situación final 
Ux = [];
Uy = [];
for i= 1:64
  for j = 1:11
    Ux(i,j)=(-sin(pi*vv(i,j)).*sin(vv(i,j)))./(20.*uu(i,j));
    Uy(i,j)=((sin(pi.*vv(i,j))).*cos(vv(i,j)))./(20.*uu(i,j));
    Uxx(i,j) = xx(i,j)+Ux(i,j);
    Uyy(i,j) = yy(i,j)+Uy(i,j);
    j = j+1;
    end
  i= i+1;
end

subplot(1,2,2);
mesh(Uxx,Uyy,0*xx);
axis([-3,3,-3,3])
view(2)

4.2 Variación del volumen sobre la superficie

Con el fin de definir la variación del volumen de la placa tras la percusión, usamos el operador divergencia. Aplicando este operador al campo vectorial [math]\vec u[/math], obtenemos la función escalar divergencia de [math]\vec u[/math]: [math]\nabla \cdot \vec u =\frac{\pi\cos(\pi \theta)}{20\rho^2}\ [/math]. Representando esta función mediante el código:

i=1;
j=1;
for i = 1:64			% Cálculo de la divergencia en cada punto
  for j = 1:11
    divergente(i,j)= (pi.*(cos(pi.*vv(i,j))./(20.*(uu(i,j).^2)));
    j= j+1;
  end
  i = i+1;
end
hold on
mesh(xx,yy,divergente); 
contour(xx,yy,divergente,20);
hold off

subplot(1,2,2);
contour(xx,yy,divergente,20); colorbar;

Se observa que el máximo aumento de volumen se produce en zonas radiales, que coinciden con los lugares hacia los que se orienta la mayoría de los vectores del campo [math]\vec u[/math]. Obviamente, al haber zonas con aumento máximo de volumen, también existen zonas cuya disminución de volumen es máxima. Estas zonas se ubican alternadas entre sí, distinguiéndose, tres áreas principales en las que se concentrará la materia, y tres en las que se reducirá la cantidad de materia.

4.3 Tendencia a la rotación del desplazamiento

Se define el rotacional como la tendencia de un campo vectorial a inducir rotación alrededor de un punto. En nuestro campo de desplazamiento [math]\vec u[/math], el rotacional será: [math]\nabla \times \vec u =\frac{-sin(\pi \theta)}{10\rho^2}\vec g_{z}\ [/math] Aplicando este rotacional a cada punto de la superficie, podremos analizar en qué puntos de la misma encontramos un mayor rotacional, es decir, dónde la percusión inducirá una mayor rotación.

rotacional = zeros(64,11);
i=1;
j=1;
for i = 1:64
  for j = 1:11
    rotacional(i,j)= (sin(pi.*vv(i,j))./(uu(i,j).^2)); %% módulo del rotacional
    j= j+1;
  end
  i = i+1;
end
subplot(1,2,1)
quiver3(xx,yy,xx*0,xx*0,xx*0,-rotacional); % representación campo vectorial

subplot(1,2,2)
mesh(xx,yy,rotacional); colorbar;  % módulo del campo rotacional

4.4 Tensiones provocadas por la perturbación

En un medio elástico lineal, isótropo y homogéneo los desplazamientos permiten escribir el tensor de tensiones [math]\sigma_{ij}[/math] a través de la fórmula: [math]\sigma_{ij}=\lambda \nabla \cdot \vec u \delta_{ij} + 2\mu \epsilon_{ij},[/math] Dónde [math]\lambda[/math] y [math]\mu[/math] son los conocidos como coeficientes de Lamé que dependen de las propiedades elásticas de cada material. En este caso, tomando [math]\lambda=\mu=1[/math], dibujaremos las tensiones normales en la dirección [math]\vec g_{\rho}[/math], es decir [math]\vec g_{\rho} \cdot \sigma \cdot \vec g_\rho[/math] y las tensiones normales en la dirección [math]\vec g_\theta/\rho[/math], es decir [math]\vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho[/math]

[math]\sigma_{ij}=\begin{bmatrix} \frac{\pi\cos(\pi \theta)}{20\rho^2} & \frac{-sin(\pi \theta)}{40}(\frac{1}{\rho}+\frac{1}{\rho^3}) & 0 \\ \frac{-sin(\pi \theta)}{40}(\frac{1}{\rho}+\frac{1}{\rho^3}) & \frac{\pi\cos(\pi \theta)}{10\rho^2} & 0 \\ 0 & 0 & \frac{\pi\cos(\pi \theta)}{20\rho^2} \end{bmatrix}\\ \vec g_{\rho} \cdot \sigma \cdot \vec g_\rho = \frac{\pi\cos(\pi\theta)}{20\rho^2}\\ \vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho =\frac{\pi\cos(\pi\theta)}{10\rho^2}[/math]

%% tensiones normales en la dirección de los ejes gro y gtheta/ro
sigmaro = (pi*cos(pi.*vv))./(20.*(uu.^2));

sigmatheta = (pi*cos(pi.*vv))./10.*(uu.^2);
%% representación del mallado y de las tensiones en plano
subplot(2,2,1);
mesh(xx,yy,sigmaro); colorbar
view(2)
subplot(2,2,3);
mesh(xx,yy,sigmatheta); colorbar
view(2)
%% representación del mallado y de las tensiones en 3D
subplot(2,2,2);
mesh(xx,yy,sigmaro);
subplot(2,2,4);
mesh(xx,yy,sigmatheta);

Como podemos ver en las imágenes, las tensiones normales en la dirección [math]\vec g_{\rho}[/math] y en las de dirección [math]\vec g_\theta/\rho[/math] son mayores en las zonas que presentan una máxima deformación consecuencia de la percusión. Aun así, es importante notar que en el caso de las tensiones normales en la dirección [math]\vec g_{\rho}[/math], aumenta su intensidad según nos aproximamos al borde interior de la placa; mientras que las definidas por el vector [math]\vec g_\theta/\rho[/math], aumentan según nos acercamos al borde externo de la misma.

4.4.1 Tensiones tangenciales respecto al plano ortogonal a [math]\vec g_\rho[/math]

Calcularemos las tensiones tangenciales al plano ortogonal a [math]\vec g_{\rho}[/math], es decir, [math]|\sigma \cdot \vec g_\theta/\rho-(\vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho) \vec g_\theta/\rho|[/math]. La representación de estas, nos indicará en qué puntos las tensiones son mayores, y dado que el estudio se realiza en un medio elástico lineal, isótropo y homogéneo, en estos puntos citados las deformaciones también serán máximas.

sigmaro = (pi*cos(pi.*vv))./(20.*(uu.^2)); 	%% tensión en dirección gro
subplot(1,3,1);
mesh(xx,yy,sigmaro); 
subplot(1,3,2); 	%% tensión tangencial respecto al plano ortogonal a gro (3D)
tenstangencialplanoro=(((sin(pi.*vv))./40).*(uu.^(-1)+uu.^(-3)));
mesh(xx,yy,tenstangencialplanoro);
subplot(1,3,3); %% tensión tangencial respecto al plano ortogonal a gro (2D)
mesh(xx,yy,tenstangencialplanoro); colorbar
view(2);

4.4.2 Tensiones tangenciales respecto al plano ortogonal a [math]\vec g_\theta/\rho[/math]

Calcularemos las tensiones tangenciales al plano ortogonal a [math]\vec g_\theta/\rho[/math], es decir, [math]|\sigma \cdot \vec g_\theta/\rho-(\vec g_\theta/\rho \cdot \sigma \cdot \vec g_\theta/\rho) \vec g_\theta/\rho|[/math]. La representación de estas, nos indicará en qué puntos las tensiones son mayores, y dado que el estudio se realiza en un medio elástico lineal, isótropo y homogéneo, en estos puntos citados las deformaciones también serán máximas.

%% tensión en dirección gtheta/ro
sigmatheta = (pi*cos(pi.*vv))./10.*(uu.^2); 
subplot(1,3,1);
mesh(xx,yy,sigmatheta); 
%% tensión tangencial respecto al plano ortogonal a gtheta/ro 3D
tenstangencialplanotheta = (((sin(pi.*vv))./(40.*uu)).*(uu.^(-1)+uu.^(-3))); 
subplot(1,3,2);
mesh(xx,yy,tenstangencialplanotheta); 
%% tensión tangencial respecto al plano ortogonal a gtheta/ro 2D
subplot(1,3,3); 
mesh(xx,yy,tenstangencialplanotheta); colorbar
view(2);

Podemos concluir que en las zonas donde se produce un mayor desplazamiento, las tensiones que presenta el sólido son mayores. Asimismo se puede observar la ortogonalidad de ambas tensiones tangenciales al plano ortogonal a [math]\vec g_{\rho}[/math] y al plano ortogonal a [math]\vec g_\theta/\rho[/math] lo que confirma que las tensiones del sólido se pueden descomponer en dos componentes [math]\vec g_\theta/\rho[/math] y [math]\vec g_{\rho}[/math].