Deformaciones de una placa plana. Grupo 1-B

Trabajo realizado por estudiantes
Título	Deformaciones de una placa plana. Grupo 1-B
Asignatura	Teoría de Campos
Curso	2020-21
Autores	José Penadés Crespí; Oliver Prada Sanchidrián; Jorge Amador Domínguez; David Pallarés Carrillo
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

Participantes

José Penadés Crespí

Oliver Prada Sanchidrián

Jorge Amador Domínguez

David Pallarés Carrillo

Consideramos una placa rectangular plana (en dimensión 2) que ocupa la región (x, y) ∈ [0, 10] × [−1, 1]. En ella vamos a suponer que tenemos definidas dos cantidades físicas: la temperatura T(x, y), que viene dada por [math]T(x,y) = x^2 + (10y)^2 [/math], y los desplazamientos [math]\vec u(x,y) [/math] producidos por la acción de una fuerza determinada.

De esta forma, si definimos [math]\vec r_{0}(x,y) [/math] el vector de posición de los puntos de la placa antes de la deformación, la posición de cada punto (x, y) de la placa después de la deformación viene dada por: [math]\vec r(x,y) = \vec r_{0}(x,y) + \vec u(x,y) [/math].

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 [math] \vec u(x,y)= xf(y) \vec i − (\frac{x^2}{200}) \vec j [/math], donde f(y) es una cierta función que no conocemos.

1 INTRODUCCIÓN

Este trabajo se centra en el estudio físico de una placa rectangular en dos dimensiones. Se ha realizado analizando temperatura y deformaciones tratándolos como campos vectoriales. Para llevar a cabo todos los cálculos necesarios se ha empleado el software científico de visualización Octave.

2 MALLADO

Para representar el mallado del sólido se ha dibujado sobre un rectángulo (x,y) ∈ [−0.5; 10.5] × [−1.5; 1.5] y como paso de muestreo h = 1/10 para las variables x e y.

% Región de la placa  
 h=0.1;
 u=0:h:10;                  
 v=-1:h:1;
 % Matriz de u y v
 [Mx,My]=meshgrid(u,v);
 %Mallado
 mesh(Mx,My,0*Mx);
 %Región de dibujo
 axis([-0.5,10.5,-1.5,1.5]);
 axis equal
 title ('Placa')
 view(2)
 xlabel('Eje X')
 ylabel('Eje Y')

3 CURVAS DE NIVEL

Representamos la variación de temperaturas a partir de la siguiente fórmula: [math]T(x,y) = x^2 + (10y)^2 [/math]. La gráfica obtenida en Octave muestra una escala de colores que varía en función de la temperatura, desde temperaturas más frías en color azul, hasta temperaturas más cálidas en colores naranjas. Por lo tanto, la temperatura será máxima en los puntos de la gráfica con los colores más cálidos, en nuestro caso, cuando la x aumenta y la y se aproxima a 1 y -1.

h=0.1;
u=0:h:10;                  
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
%Definimos función temperatura
Temp= (Mx.^2+100.*My.^2);    
%Representamos líneas de nivel
contour(Mx,My,Temp);         
axis([-0.5,10.5,-1.5,1.5]);
axis equal
xlabel('Eje X')
ylabel('Eje Y')
colorbar
title ('Distribución de la temperatura')

4 GRADIENTE DE TEMPERATURA DE LA PLACA

Una vez se ha obtenido la temperatura, procedemos a obtener como de rápido varía la temperatura cuando nos movemos desde un punto de la placa en una dirección determinada, es decir, el gradiente de temperatura. Derivamos respecto a cada variable, la función de la temperatura expresada en el punto anterior. Tras su cálculo, terminamos representando el campo vectorial sobre las curvas de nivel del campo escalar y observar su variación gráficamente. Como el gradiente muestra la dirección de máxima variación en cada punto, los vectores deben ser perpendiculares a las curvas de nivel, como se puede observar en la imagen.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
Temp= (Mx.^2+100.*My.^2);
contour(Mx,My,Temp);
%Añadimos al dibujo de las líneas de nivel del gradiente de temperatura
hold on
[Px,Py]=gradient (Temp,0.1,0.1);    
quiver(u,v, Px, Py)              
axis([-0.5,10.5,-1.5,1.5]);
axis equal
hold off
xlabel('Eje X')
ylabel('Eje Y')
title ('Gradiente de la Temperatura')

5 CAMPO DE DESPLAZAMIENTOS

Consideramos ahora el campo de desplazamientos [math] \vec u [/math] definido como: [math] \vec u(x,y)= xf(y) \vec i − (\frac{x^2}{200}) \vec j [/math] .

Para calcular f(y), se nos dice que los puntos situados en el eje y = 0 no sufren desplazamiento en la dirección [math] \vec i [/math] y que el rotacional de [math] \vec u [/math] es igual a [math] \frac{3x}{100}. [/math] Por tanto si calculamos el rotacional de u y lo igualamos a [math] \frac{3x}{100} [/math] obtendremos f(y). En este caso el resultado obtenido es que f(y) es igual a [math] \frac{y}{50} [/math] .

Volviendo a la expresión inicial y sustituyendo, ahora podemos decir que [math] \vec u= \frac{xy}{50}\vec i - \frac{x^2}{200}\vec j [/math] .

El desplazamiento se produce en la parte derecha de la placa siguiendo una dirección vertical de sentido negativo y una dirección horizontal en sentido positivo dejando la placa deformada. La parte izquierda de la placa no sufre movimiento alguno.

6 CAMPO DE VECTORES

A continuación se representa mediante un gráfico realizado en Octave, la representación del campo de vectores en los puntos del mallado del sólido.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
%Definimos vector desplazamiento  
Ux= Mx.*My/50;                
Uy= -Mx.^2 /200;
%Representamos campo vectorial
quiver(Mx,My,Ux,Uy);           
axis([-0.5,10.5,-1.5,1.5]);
axis equal
xlabel('Eje X')
ylabel('Eje Y')
view(2)
title ('Campo de desplazamientos')

7 SÓLIDO ANTES Y DESPUÉS DEL DESPLAZAMIENTO

En este apartado, procederemos a mostrar gráficamente el antes y el después del desplazamiento de la placa, junto con una superoposición de las imágenes con la función plot3 para que se pueda apreciar visualmente la deformación.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
Ux= Mx.*My/50 ;
Uy= -Mx.^2 /200 ;
figure(1)
%ANTES
subplot(2,2,1)
mesh(Mx,My,0*Mx);
axis([-0.5,10.5,-1.5,1.5]);
axis equal
title('Antes del desplazamiento');
view(2)
xlabel('Eje X')
ylabel('Eje Y') 
%DESPUÉS
subplot(2,2,2)
mesh(Ux+Mx ,Uy+My,0*Ux);
axis([-0.5,10.5,-1.5,1.5]);
axis equal 
title('Después del desplazamiento')
view(2)
xlabel('Eje X')
ylabel('Eje Y') 
%COMPARACIÓN
subplot(2,2,3)
plot3(Mx,My,Mx*0);
hold on 
plot3(Ux+Mx,Uy+My,0*Ux);
hold off
view(2)
axis([-0.5,10.5,-1.5,1.5]);
axis equal
title ('Comparación')
hold off
xlabel('Eje X')
ylabel('Eje Y')

8 DIVERGENCIA DEL CAMPO VECTORIAL SOBRE NUESTRA PLACA

Una vez observada la deformación de nuestra placa, estudiaremos desde un punto de vista muy parecido a como veíamos el campo escalar de temperatura en nuestra placa, pero con la deformación de la misma. Usaremos un sistema de colores, igual que en el campo escalar pero esta vez, del color más frío al más cálido, indicará menos deformación a más deformación, respectivamente, así como la indicación de la divergencia máxima y mínima de la placa. Por tanto, se puede intuir debido a la relación del título con esta explicación, que la divergencia, da esta información que acabo de destacar, es decir, hace referencia a cuanto se deforma la placa en cualquier punto de la misma.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
figure(2)
%Representamos la divergencia
div= My./50;           
surf(Mx,My, div);
axis equal
view(2)
xlabel('Eje X')
ylabel('Eje Y')
colorbar
title ('Divergencia')

9 ROTACIONAL DE [math]\vec u [/math]

Con la ayuda del programa Octave, se observa que, se ha calculado el capo rotacional del campo desplazamiento ?u en los puntos del sólido y se ha representado como se observa en la gráfica, qué puntos sufren mayor rotacional.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
%Definimos rotacional 
Rot= -3*Mx/100;              
%Representamos rotacional
surf(Mx,My,Rot);             
axis([-0.5,10.5,-1.5,1.5]);
axis equal
view(2)
colorbar
xlabel('Eje X')
ylabel('Eje Y')
title ('Rotacional')

10 TENSOR DE TENSIONES

Definamos �(~u) = (∇~u + ∇~ut )/2, la parte simétrica del tensor gradiente de ~u conocido como tensor de deformaciones. En un medio elástico lineal, isótropo y homogéneo los desplazamientos permiten escribir el tensor de tensiones σij a través de la fórmula: σ = λ∇ · ~u 1 + 2µ�, donde λ y µ son los conocidos como coeficientes de Lamé que dependen de las propiedades elásticas de cada material. A pesar de que los desplazamientos son planos (es decir ~u no tiene componente en la dirección de ~k) las tensiones no tienen por qué ser planas y puede haber tensiones en la dirección ortogonal al plano de la placa. Tomando λ = µ = 1, dibujar las tensiones normales en la dirección que marca el eje ~i, es decir ~i · σ ·~i, las tensiones normales en la dirección que marca el eje ~j, es decir ~j · σ · ~j y las correspondientes al eje ~k, es decir ~k · σ · ~k.

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
%tensiones normales en la dirección que marca el eje i
subplot (2,2,1)
tx=0;
ty=3.*My/50;
quiver(Mx,My,tx,ty);
xlabel('Eje X')
ylabel('Eje Y')
axis([-0.5,10.5,-1.5,1.5]);
axis equal
title ('Tensiones normales a la dirección i')
colorbar

%tensiones normales en la dirección que marca el eje j
 subplot(2,2,2)
 tx=0;
ty=My/50;
quiver(Mx,My,tx,ty);
xlabel('Eje X')
ylabel('Eje Y')
axis([-0.5,10.5,-1.5,1.5]);
axis equal
title ('Tensiones normales a la dirección j')
colorbar

%tensiones normales en la dirección que marca el eje k
 subplot (2,2,3)
tx=0;
ty=My/50;
quiver(Mx,My,tx,ty);
xlabel('Eje X')
ylabel('Eje Y')
axis([-0.5,10.5,-1.5,1.5]);
axis equal
title ('Tensiones normales a la dirección k')
colorbar

11 TENSORES TANGENCIALES

Al realizar el cálculo de las tensiones tangenciales de nuestra placa respecto al plano ortogonal a ?i, es decir |σ ·~i − (~i · σ ·~i)~i|, hemos obtenido lo siguiente:

[math] [σ \vec i - (\vec i σ \vec i)]=[σ \vec i - \frac{3y}{50}\vec i] = [\frac{3y}{50}\vec i + \frac{x}{100}\vec j - \frac{3y}{50}\vec i]= \frac{x}{100}\vec j [/math]

12 TENSIÓN DE VON MISES

La tensión de Von Mises (o el Esfuerzo) es un escalar obtenido de la combinación de los Esfuerzos Principales en un momento dado para determinar en qué puntos ocurre el esfuerzo en el eje X, Y y Z y provoca el fallo. Este método de cálculo se utiliza para medir el esfuerzo y las distribuciones de tensión dentro de un material dúctil. Para el calculo es necesario calcular los autovalores de la matriz de tensiones en cada punto del sólido. La matriz del tensor de tensiones se define como σ , puede calcularse fácilmente a partir de las tensiones principales (autovalores) del tensor tensión en un punto del sólido deformable, mediante la expresión, Von Mises = [math] \sqrt{\frac{(\sigma_1-\sigma_2)^2-(\sigma_2-\sigma_3)^2-(\sigma_3-\sigma_1)^2}{2}} [/math]

A continuación se incluye el cálculo realizado en Octave y su gráfico:

h=0.1;
u=0:h:10;
v=-1:h:1;
figure (3)
[Mx,My]=meshgrid(u,v);
f1=inline('3*My./50','Mx','My');
f12=inline('1/100.*Mx','Mx','My');
f13=inline('0.*Mx','Mx','My');
f21=inline('1/100.*Mx','Mx','My');
f2=inline('My./50','Mx','My');
f23=inline('0.*Mx','Mx','My');
f31=inline('0.*Mx','Mx','My');
f32=inline('0.*Mx','Mx','My');
f3=inline('My./50','Mx','My');
%Bucle  
for j= 1:length(u)                             
    for i= 1:length(v)
        F1=f1(Mx(i,j),My(i,j));                       
        F12=f12(Mx(i,j),My(i,j));
        F13=f13(Mx(i,j),My(i,j));
        F21=f21(Mx(i,j),My(i,j));
        F2=f2(Mx(i,j),My(i,j));                       
        F23=f23(Mx(i,j),My(i,j));
        F31=f31(Mx(i,j),My(i,j));
        F32=f32(Mx(i,j),My(i,j));
        F3=f3(Mx(i,j),My(i,j));               
        T=[F1 F12 F13;F21 F2 F23;F31 F32 F3];
        AV=eig(T);
        vmax=sqrt(((AV(1)-AV(2))^2+((AV(2)-AV(3))^2+((AV(3)-AV(1))^2))*0.5));
        Mz(i,j)=vmax;
    end
end
surf(Mx,My,Mz)
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Eje Z')
title('Tensión de Von Mises')

13 CAMPO DE FUERZAS

Para calcular el campo de fuerzas F~ que actúa sobre la placa (y que son las causantes del desplazamiento observado) se aproxima usando la ecuación de la elasticidad lineal F~ = −∇ · σ = −∂σji/∂xj.

[math] - \frac{∂σji}{∂xj}\vec e= \vec F= \frac{\frac{∂3y}{50}}{∂y}\vec i- \frac{\frac{∂x}{100}}{∂y}\vec i- \frac{∂0}{∂z}\vec i- \frac{\frac{∂x}{100}}{∂x}\vec j- \frac{\frac{∂y}{50}}{∂y}\vec j- 0 -0 -0- \frac{\frac{∂y}{50}}{∂y}\vec k = - \frac{1}{100}\vec j - \frac{1}{50}\vec j - \frac{3}{100}\vec j = F = - \frac{3}{100}\vec j [/math]

h=0.1;
u=0:h:10;
v=-1:h:1;
w=zeros(size(u));
[Mx,My,Mz]=meshgrid(u,v,w);
X=0;
Y=-3/10;
Z=0;
quiver3(Mx,My,Mz,X,Y,Z)
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Eje Z')
title('Campo de fuerzas')

14 MASA TOTAL DE LA PLACA

En el cálculo de la masa de la placa, se ha tomado la ecuación de densidad: [math] d(x, y) = 1 + xy log(1 + x + y^2). [/math]

Para ello calculamos la integral

[math]\int_{A} d(x,y) dA=\int_{0}^{10}\int_{-1}^{1} d(x,y) dydx=\int_{0}^{10}\int_{-1}^{1} d(x, y) = 1 + xy log(1 + x + y^2) [/math]

El resultado de la integral, y por tanto la masa de la placa es de 20 unidades

h=0.1;
u=0:h:10;
v=-1:h:1;
[Mx,My]=meshgrid(u,v);
%Límites de integración
a=0;                               
b=10;
c=-1;
d=1;
%Número de puntos
N1=(b-a)/h;                        
N2=(d-c)/h;
%Nuestras función
f=1+Mx.*My.*log10(1+Mx+My.*My);    
%Vectores de masa
m1=ones(N1+1,1);                   
m1(1)=0.5;
m1(N1+1)=0.5;
m2=ones(N2+1,1);
m2(1)=0.5;
m2(N2+1)=0.5;
%Resultado
MASA= h^2*m2'*f*m1 
fprintf('El valor de la masa es %.4f', MASA)