Revisión actual del 19:45 5 dic 2014

1 ENUNCIADO

Consideramos una placa plana (en dimensión 2) que ocupa la región comprendida entre las parábolas P1 : 18y-􀀀81x²2-􀀀1 = 0 y P2 : 2y+x²2-􀀀1 = 0. Para representarla usaremos un sistema de coordenadas adaptado a la geometría que nos dan:

• x=uv
• y=(1/2)(u²-v²)

con u y v definidas en (u,v) 2 [1/3,1]x[􀀀-1,1].En ella vamos a suponer que tenemos definidas dos cantidades físicas: la temperatura T(u,v), que depende de las dos coordenadas curvilíneas (u,v), y los desplazamientos u(x,y) producidos por la acción de una fuerza determinada. De esta forma, si definimos r0(u,v) el vector de posición de los puntos de la placa antes de la deformación, la posición de cada punto (u,v) de la placa después de la deformación viene dada por:

• r(u,v) = r₀(u,v) + u(u,v)

Vamos a suponer que la fuerza aplicada sobre la placa ha provocado un desplazamiento de los puntos de la misma dado por el vector de desplazamientos u(u,v) = a (b . r₀) donde a y b son vectores dados. En este trabajo supondremos lo siguiente: a=g_u/|g_u| y b=(-4).(g_u/|g_u|)

1.1 Representación gráfica de la placa

En primer lugar, se obtiene el gráfico que representa la placa que consideraremos a lo largo de todo este trabajo, que corresponde a la región comprendida entre dos parábolas dadas
(P1 : 18y-􀀀81x²-􀀀1 = 0 y P2 : 2y+x²􀀀-1 = 0).El siguiente programa realiza el mallado correspondiente a (u,v) en el intervalo [-1,1]x[-1,1], y con un paso de 1/20 (=h).

h=(1/20);
u=(1/3):h:1;
v=-1:h:1;
[uu,vv]=meshgrid(u,v);
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
figure(1)
subplot(1,2,1);
mesh(xx,yy,0*xx);
axis([-1,1,-1,1]);
view(2);
subplot(1,2,2);
mesh(xx,yy,0*xx);
axis([-1,1,-1,1]);

2 Lineas coordenadas y vectores de la base natural

Definimos las líneas coordenadas como curvas que se obtienen al fijar dos de las tres coordenadas y hacer variar la tercera. En este apartado dibujaremos dichas líneas y los vectores de la base natural,que en nuestro caso son los siguientes: [math] \vec{g_u}=v\vec{e_1} +u \vec{e_2}[/math] : [math] \vec{g_v}=u\vec{e_1} -v \vec{e_2}[/math] : [math]\vec{g}_w=\vec{e_3} [/math]

(Con el siguiente programa obtenemos dos gráficos, en el primero de ellos se aprecia g_u aplicado en cada punto de la malla, y en el segundo aparece aplicado g_v )

h=(1/20);
u=(1/3):h:1;
v=-1:h:1;
[uu,vv]=meshgrid(u,v);
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
guu=vv; guv=uu; gvu=uu; gvv=-vv;
figure(2)
subplot(1,2,1);
quiver(xx,yy,guu,guv);
axis([-1.6,1.6,-0.8,0.8]);
subplot(1,2,2);
quiver(xx,yy,gvu,gvv);
axis([-1.6,1.6,-0.8,0.8]);

Aplicando ambas componentes de la base natural sobre cada punto de la malla, todo en un mismo gráfico, se obtiene:

figure (3)
hold on
quiver(xx,yy,guu,guv);
quiver(xx,yy,gvu,gvv);
axis([-1.6,1.6,-0.8,0.8]);
hold off

INTERPRETACIÓN:Como podemos observar, los vectores g_u y g_v son claramente ortogonales entre sí en todos los puntos de la malla.

3 Campo de temperaturas

Tomamos el campo de temperaturas dado por la función T(x,y)=(8-y²+2y).e^-(x2)

3.1 Máxima temperatura del campo

h=(1/20);
u=(1/3):h:1;
v=-1:h:1;
[uu,vv]=meshgrid(u,v);
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
f=(8-(yy.^2)+(2.*yy)).*exp(-(xx.^2));
figure(102)
subplot(1,2,1);
surf(xx,yy,f);
view(2);
colorbar
axis([-1,1,-1,1]);
subplot(1,2,2);
surf(xx,yy,f)
axis([-1,1,-1,1])

INTERPRETACIÓN:En este gráfico vemos la actuación del foco de calor sobre la placa plana, observándose como varía dependiendo de la distancia a dicho foco. A la hora de obtener el máximo de temperatura, hemos calculado las derivadas parciales del campo escalar y las hemos igualado a cero. Con esto hemos observado un máximo en (0,1), punto que no entra dentro del dominio en estudio. Observando el gráfico se aprecia que el punto más cercano al mencionado anteriormente es (0,0,5).

3.2 Líneas de nivel del campo

h=(1/20);
u=(1/3):h:1;
v=-1:h:1;
[uu,vv]=meshgrid(u,v);
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
f=(8-(yy.^2)+(2.*yy)).*exp(-(xx.^2));
figure(104)
hold on
contour(xx,yy,f,20);
view(2);
axis([-1,1,-1,1]);
hold off

3.3 Variación de temperatura

Para obtener dicha variación de temperaturas, se calcula el gradiente de la función del campo escalar, que nos indica la dirección en la que dicho campo varía con mayor rapidez y que en este caso es:

[math]\nabla T(x,y)= \frac{\partial T}{\partial x} \vec{e_1}+\frac{\partial T}{\partial y}\vec{e_2}=[-2xe^{-x^2}(8-y^2+2y)]\vec{e_1}+[(-2y+2)e^{-x^2}]\vec{e_2}[/math],

h=(1/20);
u=(1/3):h:1;
v=-1:h:1;
[uu,vv]=meshgrid(u,v);
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
f=(8-(yy.^2)+(2.*yy)).*(exp(-(xx.^2)));
Tx=-(8-(yy.^2)+(2.*yy)).*((2.*xx).*exp(-(xx.^2)));
Ty=((-2.*yy+2).*exp(-(xx.^2)));
figure(152)
subplot(1,3,1),
quiver(xx,yy,Tx,Ty)  
axis([-1,1,-1,1]);
subplot(1,3,2)
hold on
quiver(xx,yy,Tx,Ty)     
contour(xx,yy,f,10)
axis([-1,1,-1,1]);
axis equal
hold off

INTERPRETACIÓN: Se puede observar, por la propia definición del gradiente, que dicho campo vectorial es perpendicular a las curvas de nivel de nuestra placa plana.

4 Campo de desplazamientos

Nos disponemos a estudiar como actúa sobre la placa el vector de desplazamientos [math]\vec{u}[/math], siendo el nuevo vector de posiciones de cada punto de la malla la suma de su posición inicial ([math]\vec{r_0)}[/math] mas dicho vector [math]\vec{u}[/math], ambos en función de las coordenadas curvilíneas (u,v).

Desarrollando el vector de las nuevas posiciones tras el desplazamiento obtenemos:

[math]\vec{u}= \frac{\vec{g}_u}{|\vec{g_u}|}(-4 \frac{\vec{g}_u}{|\vec{g_u}|} \vec{r_o})= \frac{\vec{g}_u}{|\vec{g_u}|^2}(-4 \vec{g_u} \vec{r_o})=-2u\vec{g_u}=-2uv \hat{e_1} -2u^2 \hat{e_2}[/math]

4.1 Representación campo de vectores

h=1/20;                 
u=1/3:h:1;              
v=-1:h:1;               
[uu,vv]=meshgrid(u,v);    
xx=uu.*vv;                
yy=1/2*(uu.^2-vv.^2);  
figure(76)
Ux=-2.*uu.*vv;
Uy=-2.*uu.^2; 
quiver(xx,yy,Ux,Uy)
axis equal
view(2)

4.2 Sólido antes y después del desplazamiento

u=1/3:h:1;            
v=-1:h:1;               
[uu,vv]=meshgrid(u,v);    
xx=U.*V;                 
yy=1/2*(U.^2-V.^2);  
Ux=-2.*uu.*vv;
Uy=-2.*uu.^2;
figure(5)
Dx=xx+Ux; 
Dy=yy+Uy; 
subplot(1,3,1)
mesh(xx,yy,0*xx)
axis([-1,1,-2,1]) 
view(2) 
title('Sólido inicial')
axis equal      
view(2)                
subplot(1,3,2)
mesh(Dx,Dy,0*Dx)
axis([-1,1,-2,1]) 
view(2)  
title('Sólido con el desplazamiento')
axis equal    
view(2)               
subplot(1,3,3)
hold on
mesh(Dx,Dy,0*Dx)
quiver(xx,yy,Ux,Uy)
title('Sólido con el desplazamiento y el campo del desplazamiento')
axis equal    
view(2)  
hold off

INTERPRETACIÓN:En la segunda gráfica se observa el desplazamiento que sufre la placa inicial tras aplicar el vector [math]\vec{u}[/math], en la cual se aprecia que los puntos que sufren un mayor desplazamiento son aquellos situados en los extremos inferiores. En la tercera gráfica superponemos las dos anteriores para una visión más clara de la deformación sufrida.

4.3 Variación de volumen producida por desplazamiento

Para calcular esta variación de volumen local calculamos la divergencia del campo de desplazamiento, que es el siguiente:

[math]\nabla· \vec{u}= \frac{1}{ \sqrt{g} } \frac{\partial \sqrt{g} u^{i} }{\partial u^i}= \frac{1}{u^2+v^2}(-6u^2-2v^2)= \frac{-(6u^2+2v^2)}{u^2+v^2} [/math]

Al encontrarnos con una divergencia negativa, sabemos que se trata de un sumidero, es decir, la diferencia entre el flujo saliente y entrante es menor que cero; en el caso de haber salido positiva, se trataría de una fuente.

h=(1/20);                  
u=[(1/3):h:1];               
v=[-1:h:1];            
[uu,vv]=meshgrid(u,v); 
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
div=(-6*uu.^2-2*vv.^2)./(uu.^2+vv.^2);
subplot(1,2,1)
surf(xx,yy,div)
view(2)
subplot(1,2,2)
surf(xx,yy,div);      
colorbar 
axis([-1.2,1.2,-0.6,0.6])   
axis equal

INTERPRETACIÓN:Como hemos dicho anteriormente, este gráfico representa las variaciones de volumen y se corrobora lo comentado en la gráfica del apartado anterior: los puntos de la placa que sufren una mayor cambio de volumen son aquellos situados en los extremos inferior. Por lo tanto son éstos los puntos de mayor divergencia.

4.4 Tendencia del campo a inducir rotación alrededor de un punto

Para calcular esta tendencia se calcula el rotacional del campo de desplazamientos, y lo que obtenemos en nuestro caso es:

[math]\nabla \times\vec{u}= \frac{1}{ \sqrt{g} } \begin{bmatrix} g_{u} & g_{v} & g_{w} \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\ u_{u} &u_{v}&u_{w}\end{bmatrix} = \frac{1}{ u^2+v^2 } \begin{bmatrix} g_{u} & g_{v} & g_{w} \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\ -2u(u^2+v^2) &0&0\end{bmatrix}=\frac{u^2+v^2}{ u^2+v^2 } \begin{bmatrix} g_{u} & g_{v} & g_{w} \\ \frac{\partial}{\partial u} & \frac{\partial}{\partial v} & \frac{\partial}{\partial w} \\-2u&0&0\end{bmatrix}=0[/math]

Obtenemos que el rotacional es nulo (irrotacional), lo cual implica que la divergencia es distinta de cero, como se ha comprobado anteriormente.

4.5 Tensor de tensiones

Se define el tensor de tensiones como [math]σ=λ\nabla·\vec{u}1+2μЄ[/math] donde λ y μ son los denominados coeficiente de Lamé (dos constantes elásticas que caracterizan por completo el comportamiento elástico de un sólido isótropo para pequeñas deformaciones), que en este caso los hemos considerado iguales a la unidad.

Para comenzar a calcular el tensor de tensiones empezamos por obtener: [math]\epsilon (\vec{u})=\frac{ \nabla \vec{u}+ \nabla \vec{u}^t}{2}[/math] para lo que debemos calcular tanto [math]\nabla\vec{u}[/math] como [math] \nabla \vec{u}^t[/math].

[math]\nabla\vec{u}=\begin{bmatrix}-2- \frac{2u^2}{u^2+v^2} & \frac{-2uv}{u^2+v^2} & 0 \\\frac{2uv}{u^2+v^2} & \frac{-2u^2}{u^2+v^2} & 0\\0&0&0 \end{bmatrix}[/math] [math] \nabla \vec{u}^t=\begin{bmatrix}-2- \frac{2u^2}{u^2+v^2} & \frac{2uv}{u^2+v^2} & 0 \\\frac{-2uv}{u^2+v^2} & \frac{-2u^2}{u^2+v^2} & 0\\0&0&0 \end{bmatrix} [/math]

y así obtenemos

[math] Є(\vec{u})=\begin{bmatrix}\frac{-8u^2-4v^2}{2(u^2+v^2)} & 0 & 0 \\\ 0 & \frac{-4u^2}{2(u^2+v^2)} & 0\\0&0&0 \end{bmatrix} [/math]

Para finalizar, con la fórmula anteriormente definida y habiendo calculado la divergencia, calculamos el tensor de tensiones.

[math] \sigma ^ij=\frac{-6u^2-2v^2}{u^2+v^2} \begin{bmatrix}1&0&0 \\0&1&0\\0&0&1 \end{bmatrix} +2 \begin{bmatrix} \frac{-8u^2-4v^2}{2(u^2+v^2)}&0&0 \\0&\frac{-8u^2-4v^2}{2(u^2+v^2)}&0\\0&0&0 \end{bmatrix}=\begin{bmatrix}\frac{-14u^2-6v^2}{u^2+v^2} & 0 & 0 \\\ 0 & \frac{-10u^2-2v^2}{u^2+v^2} & 0\\0 &0& \frac{-6u^2-2v^2}{u^2+v^2} \end{bmatrix} [/math]

A continuación dibujamos las tensiones normales en las direciones [math] \vec{g}_u [/math] y [math] \vec{g}_v [/math].

4.5.1 Tensiones normales en la dirección de g_u

h=(1/20);                  
u=[(1/3):h:1];               
v=[-1:h:1];            
[uu,vv]=meshgrid(u,v); 
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
figure(90)
Eu=(-14*uu.^2-6*vv.^2)./(uu.^2+vv.^2); 
figure(90)
subplot(1,2,1);
surf(xx,yy,Eu)      
axis([-1,1,-1,1])      
view(2)
subplot(1,2,2)
surf(xx,yy,Eu)      
axis([-1,1,-1,1])

4.5.2 Tensiones normales en la dirección de g_v

h=(1/20);                  
u=[(1/3):h:1];               
v=[-1:h:1];            
[uu,vv]=meshgrid(u,v); 
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
Ev=(-10*uu.^2-2*vv.^2)./(uu.^2+vv.^2); 
figure(92)
subplot(1,2,1);
surf(xx,yy,Ev)      
axis([-1,1,-1,1])      
view(2)
subplot(1,2,2)
surf(xx,yy,Ev)      
axis([-1,1,-1,1])

4.6 Tensión de Von Mises

La tensión de Von Mises es un escalar proporcional a la energía de deformación elástica de distorsión que puede expresarse en función de las componentes del tensor tensión, en particular admite una expresión simple en función de las tensiones principales. Se define con la fórmula

[math] \sigma _{VM}= \sqrt{ \frac{( \sigma _1- \sigma _2)^2+(\sigma_2-\sigma_3)^2+(\sigma_3-\sigma_1)^2 }{2} }[/math],

donde [math] \sigma_1 \sigma_2 [/math] y [math] \sigma_3 [/math] son los autovalores del tensor σ (tensiones principales).

h=(1/20);                  
u=[(1/3):h:1];               
v=[-1:h:1];            
[uu,vv]=meshgrid(u,v); 
xx=uu.*vv;
yy=(1/2).*((uu.^2)-(vv.^2));
figure(95)
e1=(-14*uu.^2-6*vv.^2)./(vv.^2+uu.^2); 
e2=(-10*uu.^2-2*vv.^2)./(vv.^2+uu.^2);   
e3=(-6*uu.^2-2*vv.^2)./(vv.^2+uu.^2);    
VonMises=sqrt(((e1-e2).^2+(e2-e3).^2+(e3-e1).^2)./2);
subplot(1,2,1);  
surf(xx,yy,VonMises) ;       
axis([-1,1,-1,1]);      
axis equal;
view(2);
subplot(1,2,2);
surf(xx,yy,VonMises) ;       
axis([-1,1,-1,1]);      
axis equal;

INTERPRETACIÓN:Como observamos en la imagen, la zona más propensa a deformarse plásticamente es la del eje de simetría central de la placa, encontrándose en esa zona los puntos que sufren una mayor tensión de Von Mises.

5 Masa de la placa

La densidad de la placa se define en el enunciado por [math] d (x,y) = xy e^{-1/x^2}[/math].

Para calcularla, ya que es imposible hacerlo a mano, utilizaremos un programa de matlab en el que aplicamos explícitamente la definición de integrales realizando una división de eje en 150 partes, que ya nos proporciona un resultado suficientemente aceptable.

h=1/150;                 
u=1/3:h:1;             
v=-1:h:1;               
[uu,vv]=meshgrid(u,v);    
X=uu.*vv;               
Y=1/2*(uu.^2-vv.^2);  
D=X.*Y.*exp(-1./X.^2); 
a=(h.^2).*D;         
MASA=sum(sum(a))

El resultado del programa es de una MASA= 9.8817e-19.