Visualización de campos escalares y vectoriales en elasticidad (Grupo 6C)
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Visualización de campos escalares y vectoriales en elasticidad. Grupo 6-C |
| Asignatura | Teoría de Campos |
| Curso | 2013-14 |
| Autores | Alejandra García-Page Acevedo, Nerea González Rivas, Jose Luis Peñaranda Ezpondaburu, Juan Antonio Rebollo Parada, Sergio Sanchez García, Clara Valle Gómez |
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
Contenido
- 1 Estado inicial del mallado
- 2 Líneas coordenadas y vectores de la base natural
- 3 Temperatura del sólido: función y gradiente
- 4 Campo de desplazamientos: visualización y movimiento
- 5 Operadores divergencia y rotacional del campo de desplazamientos
- 6 Comparación de las tensiones normales
- 7 Tensión de Von Mises
- 8 Masa de la placa
1 Estado inicial del mallado
Nos encontramos ante una placa plana bidimensional, la cual ocupa la región comprendida entre dos parábolas [math]P1=18·y−81^2−1=0[/math] y [math]P2=2·y+x^2−1=0[/math]. Esta placa ha sido sometida, en un caso de estudio, a un incremento de temperatura, mientras que en otro se ha sometido a una deformación plana. Con objeto de facilitar su cálculo, y su representación gráfica, se han usado un sistema de coordenadas curvilíneo adaptado a los límites geométricos:
[math]x=u·v[/math]
[math]y=\frac{(u^2−v^2)}{2}[/math] con u,v definidas en (u,v) ∈ [1/3,1] × [−1,1].
La representación de los puntos que conforman la placa con un paso de h:
% Mallado
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
% P1=18.*y-81.*(x.^2)-1 -----------------------------------------
% P2=2.*y+(x.^2)-1 ----------------------------------------------
% Apartado 1, dibujo --------------------------------------------
mesh(x,y,0*x) % Mallado
axis([-1,1,-1,1]) % Selecciona la región a dibujar
view(2) % Ver imagen en planta
2 Líneas coordenadas y vectores de la base natural
En nuestra placa, las lineas coordenadas y los vectores de la base natural varían según el punto de la placa que mires. Ésto se debe a que la base natural en las coordenadas cartesianas no es constante.
El programa que nos dibuja ésto es el siguiente:
%Apartado 1 - Mallado
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
% P1=18.*y-81.*(x.^2)-1;%----------------------------------------------------------------------------
% P2=2.*y+(x.^2)-1;%---------------------------------------------------------------------------------
% Apartado 2, dibujo --------------------------------------------------------------------------------
hold on
for i=1:40
a=u.*i;
b=(1/2).*((u.^2)-(i.^2));
c=v.*i;
d=(1/2).*((i.^2)-(v.^2));
end
hold on % Inicio superposición de gráficos
mesh(x,y,0*x) % Mallado completo
quiver(x,y,v,u); % Representación del primer vector de la base natural en cada punto
quiver(x,y,u,-v); % Representación del segundo vector de la base natural en cada punto
axis([-1.25,1.25,-1.25,1.25]) % Selecciona la regíon a dibujar
view(2) % Ver imagen en planta (dos dimensiones)
hold off
3 Temperatura del sólido: función y gradiente
La placa ha sido sometida a un incremento de temperatura, el cual se rige según una función escalar [math]exp(-y)[/math]. Para obtener el valor de la temperatura en cada punto, basta con aplicar la función a la malla de la placa.
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
f=exp(-(y)); % Función de Temperatura
% Apartado 3, dibujo ----------------------------------------------
subplot(1,2,1) % Dibujo de las lineas de Temperatura
contour(x,y,f,20)
axis([-1,1,-1,1])
view(2)
subplot(1,2,2) % Dibujo de la Temperatura
surf(x,y,f)
colorbar
axis([-1,1,-1,1])
view(2)
ValorMaximo=max(max(f)) % Valor máximo de la Temperatura
Asimismo, se ha buscado el incremento diferencial de la temperatura en cada punto, para lo cual se ha calculado el vector gradiente de la función temperatura.
f=exp(-y); % Función temperatura
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2));
fu =(-u.*exp(-(1/2*(u.^2-v.^2))));
fv =(v.*exp(-(1/2*(u.^2-v.^2))));
fX=v.*fu+u.*fv;
fY=u.*fu-v.*fv;
% Apartado 4, dibujo -----------------------------------
quiver(x,y,fX,fY)
axis([-1.1,1.1,-1.1,1.1])
view(2)
hold on
f=exp(-((1/2).*((u.^2)-(v.^2))));
contour(x,y,f)
axis([-1.1,1.1,-1.1,1.1])
hold off
4 Campo de desplazamientos: visualización y movimiento
La placa ha sido sometida a un campo de desplazamientos, el cual la deforma por completo. El campo viene definido por:
Tras calcular analíticamente el desarrollo de las componentes, se obtiene el resultado genérico del campo, el cual se aplica a cada punto de la placa para deformarla y obtener la placa deformada:
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices. TAMAÑO 41 14
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2));
w=u.*v.*(u.^3-2*u.*v-4*u.^2+4*v.^2);
wx=v.*(4.*v.^2).*sqrt(u.^2+v.^2); % Componente i del campo de desplazamientos
wy=u.*(4.*v.^2).*sqrt(u.^2+v.^2); % Componente j del campo de desplazamientos
% Apartado 5, dibujo ------------------------------------------------------------
subplot(1,3,2)
quiver(x,y,wx,wy)
title('Campo de deformaciones')
axis([-1.2,1.2,-1,1])
% Apartado 6, dibujo ------------------------------------------------------------
fx=x+wx;
fy=y+wy;
subplot(1,3,3)
surf(fx,fy,0*fx)
title('Placa deformada')
view(2)
subplot(1,3,1)
surf(x,y,0*fx)
title('Placa original')
view(2)
5 Operadores divergencia y rotacional del campo de desplazamientos
Para obtener información adicional acerca de la deformación de la placa, se han calculado los operadores divergencia y rotacional de la placa. El cálculo se ha realizado de forma analítica, tras lo cual se ha aplicado a la malla virtual de Matlab.
5.1 Divergencia
Con el operador genérico de la divergencia aplicado en cada punto, se ha obtenido una malla con las diferencias de flujos en cada punto de la placa.
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
% Apartado 7, dibujo ------------------------------------
div=12.*u.*v.^2./sqrt(u.^2+v.^2);
subplot(1,2,1)
surf(x,y,div);colorbar;
view(2)
subplot(1,2,2)
surf(x,y,div);colorbar;
5.2 Rotacional
Mediante el cálculo del campo del operador rotacional, se ha obtenido una visión gráfica de la tendencia de cada punto a rotar los puntos colindantes. Los valores son más altos según se aproximan a la parábola inferior de la placa, lo cual coincide con los puntos donde el campo de deformaciones es mayor.
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
% Apartado 8, dibujo ------------------------------------
rot=-8.*v.^2./(((u.^2)+(v.^2)).^(3/2));
subplot(1,2,1)
surf(x,y,rot);colorbar;
view(2)
subplot(1,2,2)
surf(x,y,rot);colorbar;
6 Comparación de las tensiones normales
Calculados la divergencia y el rotacional de [math] \vec{u} [/math], pasamos a definir el tensor de tensiones:
[math]σ=λ\nabla·\vec{u}1+2μЄ[/math] Donde [math] \epsilon [/math] es la parte simétrica de la divergencia:
[math]\epsilon (\vec{u})=\frac{ \nabla \vec{u}+ \nabla \vec{u}^t}{2}[/math] Partiendo de la hipótesis de medio elástico lineal, isótropo y homogéneo; y suponiendo los coeficientes de Lamé [math] λ=µ=1 [/math], calculamos las tensiones normales:
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2]
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2)); % Parametrización Y
Tu=28*u.*v.^2./sqrt((u.^2)+(v.^2)); %componente u del tensor de tensiones
Tv=20*u.*v.^2./sqrt((u.^2)+(v.^2)); %componente v del tensor de tensiones
% Apartado 9, dibujo -----------------------------------------------------
subplot(1,2,1)
quiver(x,y,Tu*0,Tv)
subplot(1,2,2)
quiver(x,y,Tu,Tv*0)
Se puede apreciar que tanto las tensiones en la dirección [math] \vec{g_u} [/math] como las tensiones en la dirección [math] \vec{g_v} [/math], son similares a las de la divergencia.
7 Tensión de Von Mises
La Tensión de Von Mises es una magnitud física utilizada para cuantificar la energía de distorsión. A partir de los autovalores del tensor de tensiones, es decir, las tensiones principales, se puede calcular la Tensión de Von Mises.
[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 de tensiones, [math] \sigma [/math].
h=1/20; % Paso de muestreo
uu=1/3:h:1; % Intervalo [1,2].
vv=-1:h:1; % Intervalo [0,2*pi]
[u,v]=meshgrid(uu,vv); % Matrices. TAMAÑO 41 14
x=u.*v ; % Parametrización X
y=(1/2).*((u.^2)-(v.^2));
w=u.*v.*(u.^3-2*u.*v-4*u.^2+4*v.^2);
% Apartado 10, dibujo --------------------------------------
z=[-2*v.^2,-2*u.*v];
%grad(U)=0;
a1=(28*u.*v.^2)./(((u.^2)+(v.^2)).^(3/2));
a2=(20*u.*v.^2)./(((u.^2)+(v.^2)).^(3/2));
a3=(12*u.*v.^2)./(((u.^2)+(v.^2)).^(3/2));
tvm=sqrt(((a1-a2).^2+(a2-a3).^2+(a3-a1).^2)./2);
subplot(1,2,1)
surf(x,y,tvm)
axis([-1,1,-1,1])
view(2)
subplot(1,2,2)
surf(x,y,tvm)
8 Masa de la placa
Para finalizar, se ha hecho un cálculo de la masa de la placa. Para obtenerla, se ha realizado un cálculo integral, iterado mediante el método del trapecio, lo cual devuelve un resultado de 8.0290e-005. Para obtener un resultado más preciso, se ha aumentado el número de paso hasta 1/100.
h=1/100; % Paso de muestreo.
uu=1/3:h:1; % Intervalo [1,2].
vv=-1:h:1; % Intervalo [0,2*pi].
[u,v]=meshgrid(uu,vv); % Matrices.
x=u.*v ; % Parametrización X.
y=(1/2).*((u.^2)-(v.^2));
%apartado 11
f=abs(abs(x).*exp(1).^(-1./(y.^2))); %distribución de la densidad (en valor absoluto)
a=h^2.*f;
sol=sum(sum(a))
