Campo de Temperaturas y Deformaciones para una placa 2D (Grupo 20-A)
| Trabajo realizado por estudiantes | |
|---|---|
| Título | Campo de Temperaturas y Deformaciones para una placa 2D (Grupo 20-A) |
| Asignatura | Teoría de Campos |
| Curso | |
| Autores |
|
| Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura | |
Un Campo de Temperaturas es un Campo Escalar aplicado a un objeto físico con el objetivo de representar las temperaturas propias de los puntos del objeto.
La representación visual de un campo escalar puede hacerse mediante programas como Octave y MATLAB, los cuales utilizan operaciones escalares y elementos matriciales para realizar cálculos complejos.
Contenido
- 1 Creación de la Superficie
- 2 Representación de Campos Escalares
- 3 Cálculo de ▽T y su representación
- 4 Representación del campo de vectores en los puntos del mallado del sólido
- 5 Representación de deformaciones de la placa
- 6 DIVERGENCIA DEL CAMPO DE DESPLAZAMIENTOS
- 7 Rotacional de un campo de desplazamientos
- 8 Tensiones normales
- 9 TENSIONES TANGENCIALES RESPECTO AL PLANO ORTOGONAL A [math]\ \vec i [/math]
- 10 TENSIÓN DE VON MISES
- 11 CAMPO DE FUERZAS [math]\vec{F}[/math] SOBRE LA PLACA
- 12 Módulo de desplazamiento vertical en el punto (1/2,1)
1 Creación de la Superficie
Para representar el campo de temperaturas de un sólido en MATLAB, primero ha de dibujarse un mallado que represente el objeto. Para ello se genera una malla de puntos mediante el comando meshgrid.
Imaginemos un mallado que represente los puntos interiores del sólido 2D rectangular con paso de muestreo [math](h) = 2/10[/math], ubicado en el plano cartesiano, tal que:
El código de Matlab para ello sería:
h = 0.2; %Paso de muestreo
x = -0.5:h:10.5;
y = -0.5:h:2.5;
[X,Y]= meshgrid(x,y); %Mallado
mesh(X,Y,0*X)
axis ([-1,11,-1,3])
view(2) %Vista en planta2 Representación de Campos Escalares
El campo escalar de temperaturas de un objeto físico 2D se representa al aplicar una función a cada vértice del mallado previamente generado, el cual representa el objeto. Para ello se genera una función que asigne a cada punto [math](x,y)[/math] un valor de temperatura. Siguiendo con el ejemplo anterior, una posible función de temperaturas sería:
Dicha función habría que incluirla en el código Matlab, tal que quede:
h = 0.2; %Paso de muestreo
x = -0.5:h:10.5;
y = -0.5:h:2.5;
[X,Y]= meshgrid(x,y); %Mallado
%Parametrizamos en cartesianas
T=((X-6).^2+(10*(Y-3/2)).^2); %Campo escalar de temperatura
figure
hold on
contour(X,Y,T,40) %Líneas de nivel
colorbar %Para mostrar las escalas
axis ([-1,11,-1,3])
3 Cálculo de ▽T y su representación
La divergencia de un campo escalar se expresa como la derivada parcial respecto de cada coordenada reflejado sobre la base ortonormal orientada positivamente [math]\vec{i},\vec{j},\vec{k}[/math]. Es decir:
Una vez calculado, se procede a la parametrización del gradiente en Matlab y su representación:
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
z=zeros(size(Mx)); %la z es cero por ser una placa sin espesor
mesh(Mx,My,z) %representamos de tal forma que solo se vean las líneas
t=((Mx-6).^2+(10*(My-3/2)).^2); %Campo escalar de temperatura
pcolor(Mx,My,t)
shading flat %para quitar las líneas de cuadrícula y que no se superpongan con las de nivel
contour(Mx,My,t,30,'k') %el color de las líneas de nivel sale en negro
hold on
u=(2.*Mx-12);
v=(200.*My-300);
w=(0.*Mx);
quiver3(Mx,My,z,u,v,w) %dibujamos el campo vectorial gradiente
axis([-0.5,10.5,-0.5,2.5]) %modificamos los ejes
axis equal
view(2)
hold off
El gráfico resultante representa, con líneas de nivel, las temperaturas de la placa, mientras que los vectores (azul) muestran el gradiente para una serie de puntos de referencia.
A simple vista puede parecer difícil corroborar la ortogonalidad de los vectores del gradiente con las líneas de nivel, debido a la imposibilidad de representar curvas, por lo que las líneas de nivel se representan mediante una serie de rectas cortas de referencia. Por ende, se adjunta a continuación una foto para demostrar dicha ortogonalidad.
4 Representación del campo de vectores en los puntos del mallado del sólido
El campo de vectores del mallado puede representarse mediante el uso de Octave/Matlab.
Un ejemplo del código a seguir sería el siguiente:
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
ux= 2/5.*sin(Mx);
uy=0.*My;
mesh(Mx,My,0*Mx);
hold on
quiver(Mx,My,ux,uy);
axis equal
axis([-0.5,10.5,-0.5,2.5]);
view(2)
hold off
xlabel('Eje X')
5 Representación de deformaciones de la placa
De igual forma, la representación de las deformaciones de una placa 2D siguiendo un campo vectorial [math] \vec u [/math] puede representarse siguiendo la definición de campo vectorial. Siguiendo con el ejemplo anterior, dichas deformaciones pueden representarse en programas como MatLab, ejemplificado con el código:
clc
clear
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
ux=2/5.*sin(Mx);
uy=0.*My;
Mz=zeros(size(Mx));
mesh(Mx,My,Mz);
subplot (1,2,1) %Placa antes del desplazamiento
pcolor(Mx,My,Mz)
axis equal
axis([-0.5,10.5,-0.5,2.5])
title('Antes del desplazamiento')
subplot(1,2,2) %Placa después del desplazamiento
Mxx=Mx+ux;
Myy=My+uy;
pcolor(Mxx,Myy,Mz); %Definimos la matriz post-desplazamiento
axis equal
axis([-0.5,10.5,-0.5,2.5])
title('Despues del desplazamiento')
6 DIVERGENCIA DEL CAMPO DE DESPLAZAMIENTOS
La divergencia del campo de desplazamientos también puede calcularse siguiendo el mismo procedimiento, y mediante este podemos también calcular los puntos para los cuales la divergencia será máxima, mínima y nula. Para el cálculo de la divergencia en coordenadas cartesianas, se procede a derivar el campo de desplazamientos .
Para ejemplificarlo, usamos el campo de desplazamientos hallado anteriormente:
[math]\ \vec u (x, y, z)= \frac{2}{5} \sin (x) \vec i [/math]
Por lo tanto:
[math]\ \nabla \cdot \vec u =\frac{\partial \vec u_1}{\partial x} + \frac{\partial \vec u_2}{\partial y} + \frac{\partial \vec u_3}{\partial z} = \frac{\partial (\frac{2}{5}\sin (x))}{\partial x} + \frac{\partial 0}{\partial y} + \frac{\partial 0}{\partial z} = \frac{2}{5} \cos(x) [/math]
Con esto obtenemos que la divergencia es un campo escalar que depende únicamente de la coordenada [math]\ x [/math] de cada punto
clear
clc
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
z=zeros(size(Mx)); %la z es cero por ser una placa sin espesor
mesh(Mx,My,z) %representamos de tal forma que solo se vean las líneas
g=(2/5)*cos(Mx); %Pintamos las curvas de nivel correspondientes a la temperatura
pcolor(Mx,My,g)
shading flat %para quitar las líneas de cuadrícula y que no se superpongan con las de nivel
contour(Mx,My,t,30,'k') %el color de las líneas de nivel sale en negro
colorbar
hold on
quiver3(Mx,My,z,u,v,w) %dibujamos
axis([-0.5,10.5,-0.5,2.5])
axis equal %modificamos los ejes
view(2)
hold off
La divergencia, al ser un campo escalar, puede ser representado mediante un mapa de colores 2D sobre la placa.
7 Rotacional de un campo de desplazamientos
Sea |∇ × \vec{u}| el rotacional de un campo de desplazamientos [math] \vec u = (ux, uy, uz)[/math] expresado en un sistema de coordenadas de vectores unitarios ortogonales [math]\vec{a}, \vec{b}, \vec{c}[/math], aplicado a una superficie bidimensional expresado en dicho sistema de coordenadas, se cumple que:
[math]\ \bigtriangledown \times \vec{u} = \begin{vmatrix}\vec{u} & \vec{v} &\vec{w} \\ d/da & d/db & d/dc\\ ux & uy & uz \end{vmatrix}[/math]
Por ello, para calcular el rotacional de un campo de desplazamientos, se aplica la fórmula del producto vectorial en los ejes del sistema de coordenadas.
Para el caso del sistema de coordenadas cartesiano, con ejes [math][ x, y, z ][/math], y vectores respectivos [math]\vec{i}, \vec{j}, \vec{k}[/math], se procede a calcular el rotacional.
Usando el campo [math]\vec{u} = (\vec{ux}, \vec{uy}, \vec{uz}) = (2/5sin(x), 0, 0)[/math] previo, procedemos a hacer los cálculos:
[math]\ \bigtriangledown \times \vec{u} = \begin{vmatrix}\vec{i} & \vec{j} &\vec{k} \\ d/dx & d/dy & d/dz\\ ux & uy & uz\end{vmatrix}= \begin{vmatrix}\vec{i} & \vec{j} & \vec{k} \\ d/dx & d/dy & d/dz \\ 2/5sin(x) & 0 & 0\end{vmatrix}=(0-0)\vec{i}+(0-0)\vec{j}+(0-0)\vec{k}=0 [/math]
Por lo tanto, se trata de un campo conservativo es decir, el rotacional en todos los puntos de la placa es nulo, lo que implica:
[math]\ \left | \bigtriangledown \times \vec{u} \right | = 0 [/math]
8 Tensiones normales
Sea [math]\sigma=λ\bigtriangledown \cdot \vec{u}1 + 2µe[/math] el campo de tensiones en un medio elástico, isótropo y homogéneo donde [math]e=\frac{\bigtriangledown\vec{u}+\bigtriangledown\vec{u}^t}{2}[/math] la parte simétrica del gradiente del campo vectorial [math]\vec u[/math], los coeficientes de Lamé [math]λ,µ =1[/math] y [math]\bigtriangledown \cdot \vec{u}[/math] la divergencia del campo vectorial [math]\vec u[/math] se puede proceder a los cálculos de las deformacciones en las direcciones de los vectores de la base ortonormal orientada positivamente [math]{\vec i,\vec j,\vec k}[/math]. Para ello, se acometen unos cálculos intermedios que facilitan la operación:
[math]\bigtriangledown\vec{u}=\begin{pmatrix}\frac{\partial u_{1}}{\partial x}&\frac{\partial u_{1}}{\partial y}&\frac{\partial u_{1}}{\partial z}\\\frac{\partial u_{2}}{\partial x}&\frac{\partial u_{2}}{\partial y}&\frac{\partial u_{2}}{\partial z}\\\frac{\partial u_{3}}{\partial x}&\frac{\partial u_{3}}{\partial y}&\frac{\partial u_{3}}{\partial z} \end{pmatrix}=\begin{pmatrix}\frac{2}{5}cos(x)&0&0\\0&0&0\\0&0&0\end{pmatrix}= \frac{2}{5}cos(x) \vec i \otimes \vec i[/math]
Como se observa, [math]\bigtriangledown\vec{u}[/math] es una matriz simétrica, lo que quiere decir que ella y su matriz traspuesta son idénticas, por lo tanto se representan de la misma forma. Entonces:
[math]e=\frac{\bigtriangledown\vec{u}+\bigtriangledown\vec{u}^t}{2}=\frac{\frac{2}{5}cos(x) \vec i \otimes \vec i + \frac{2}{5}cos(x) \vec i \otimes \vec i}{2}=\frac{2}{5}cos(x) \vec i \otimes \vec i[/math]
Ahora se procede al cálculo del tensor de deformaciones:
[math]\sigma=\frac{2}{5}cos(x)(\vec i \otimes \vec i + \vec j \otimes \vec j + \vec k \otimes \vec k) + 2\cdot\frac{2}{5}cos(x)(\vec i \otimes \vec i) = \frac{2}{5}cos(x)(\vec i \otimes \vec i) + \frac{2}{5}cos(x)(\vec j \otimes \vec j) + \frac{2}{5}cos(x)(\vec k\otimes \vec k) + \frac{4}{5}cos(x) (\vec i \otimes \vec i) = \frac{6}{5}cos(x) (\vec i \otimes \vec i) +\frac{2}{5} cos(x) (\vec j \otimes \vec j) +\frac{2}{5}cos(x) (\vec k \otimes \vec k) [/math]
[math]\sigma_{ij} = \begin{pmatrix} \frac{6}{5}cos(x)\ & 0 & 0 \\ 0 & \frac{2}{5}cos(x)\ & 0 \\ 0 & 0 & \frac{2}{5}cos(x)\ \end{pmatrix}[/math]
Por último se procede al cálculo de las tensiones normales a la dirección que marca cada vector de la base ortonormal orientada positivamente adoptada. Dichas tensiones son las que se encuentran en los elementos de la diagonal de la matriz [math]\sigma_{ij}[/math]:
En la dirección de [math]\vec{i}[/math]: [math]\:\:\vec i \cdot \sigma \cdot \vec i[/math] = [math]\frac{6}{5}\cos x[/math]
En la dirección de [math]\vec{j}[/math]: [math]\:\:\vec j \cdot \sigma \cdot \vec j[/math] = [math]\frac{2}{5}\cos x[/math]
En la dirección de [math]\vec{k}[/math]: [math]\:\:\vec k \cdot \sigma \cdot \vec k[/math] = [math]\frac{2}{5}cos x[/math]
Tras los cálculos, puede representarse gráficamente en Matlab. Siguiendo con el ejemplo del apartado anterior, lo hacemos utilizando el código:
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
Ti=(6/5).*cos(Mx); % Tension normal en i
Tj=(2/5).*cos(Mx); % Tension normal en j
Tk=(2/5).*cos(Mx); % Tension normal en k
figure
subplot(1,3,1)
pcolor(Mx,My,Ti)
colorbar
shading flat
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Eje Z')
axis equal
title('Tensiones normales en i')
subplot(1,3,2)
pcolor(Mx,My,Tj)
colorbar
shading flat
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Eje Z')
axis equal
title('Tensiones normales en j')
subplot(1,3,3)
pcolor(Mx,My,Tk)
colorbar
shading flat
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Eje Z')
axis equal
title('Tensiones normales en k')
9 TENSIONES TANGENCIALES RESPECTO AL PLANO ORTOGONAL A [math]\ \vec i [/math]
Sea [math]\vec{T}[/math] la tensión tangencial respecto al plano ortogonal a [math]\vec{i}[/math], que viene dado por la siguiente expresión:[math]\ \left| \sigma ·\vec i -(\vec i·\sigma·\vec i)\vec i \right| [/math]
Particularizando para nuestros datos obtenemos: [math]\frac{6}{5}\cos (x)\vec{i}-\frac{6}{5}\cos (x)\vec{i}=\vec{0}[/math]
Por lo tanto la tensión tangencial respecto al plano que contiene a los vectores [math]\vec{j}[/math] y [math]\vec{k}[/math] es nula.
10 TENSIÓN DE VON MISES
La Tensión de Von Mises es una magnitud física escalar usada en campos como la ingeniería estructural, calculada a partir de los valores de las tensiones principales de cada punto del espacio, la cual indica la tensión a aplicar a cada punto de un material para que éste inicie su comportamiento plástico.
Sea la Tensión de Von Mises [math] \sigma_{VM} [/math], calculada para un punto P de un sólido deformable, y sean las tensiones principales del tensor tensión para dicho punto [math]\sigma_1, \sigma_2, \sigma_3[/math], correspondientes a los autovalores de [math] \sigma_{ij} [/math], entonces se comprueba que la Tensión de Von Mises viene dada por la siguiente expresión:
[math] \sigma_{VM} = \sqrt{ \frac { ( \sigma_1 - \sigma_2 )^2 + ( \sigma_2 - \sigma_3 )^2 + ( \sigma_3 - \sigma_1 )^2 }{2} }[/math]
clear
clc
h=0.2;
x=-0.5:h:10.5; %Discretizamos valores de x
y=-0.5:h:2.5; %Discretizamos valores de y
[Mx,My]=meshgrid(x,y); %Creamos el mallado
Ti=(6/5).*cos(Mx); % Tension normal en i
Tj=(2/5).*cos(Mx); % Tension normal en j
Tk=(2/5).*cos(Mx); % Tension normal en k
%Autovalores de la tensión de Von Mises
Mz = ones(size(y,2),size(x,2));
T = zeros(3);
for m=1:length(y)
for n=1:length(x)
%Valores para cada componente del mallado para la función de tensión de Von Mises
a=Ti(m,n);
b=Tj(m,n);
c=Tk(m,n);
T=[a 0 0; 0,b,0; 0,0,c];
%Obtención de los autovalores de la tensión de Von Mises
[P,D]=eig(T);
Autoval1=D(1,1);
Autoval2=D(2,2);
Autoval3=D(3,3);
%Cálculo de Von Mises para cada componente
VonMises=sqrt( ((Autoval1-Autoval2)^2+(Autoval2-Autoval3)^2+(Autoval3-Autoval1)^2) * 1/2 );
Mz(m,n)=VonMises;
end
end
%Representamos Von Mises en la placa bidimensional
subplot(1,2,1)
surf(Mx,My,Mz)
colorbar
axis equal
xlabel('Eje X')
ylabel('Eje Y')
title('Tensión de Von Mises')
view(2)
%Representamos las tensiones en un tercer eje
subplot(1,2,2)
surf(Mx,My,Mz)
colorbar
axis equal
xlabel('Eje X')
ylabel('Eje Y')
zlabel('Tensión de Mises')
title('Tensión de Von Mises representada en el eje Z')
11 CAMPO DE FUERZAS [math]\vec{F}[/math] SOBRE LA PLACA
En este apartado y en adelante, supondremos que el timepo [math]t[/math] en la ecuación de la onda longitudinal es distinto de 0. Por ende, es necesario volver a calcular [math]\epsilon[/math] asumiento que [math]t[/math] no es cero. Se parte de la ecuación de elasticidad lineal que se muestra a continuación para aproximar el campo de fuerzas que actúa sobre la placa.
[math]\vec{F} = \frac{d^2\vec{u}}{dt^2} - \bigtriangledown \cdot \sigma[/math]
Donde [math]\vec{u}[/math] = [math]\frac{2}{5}\vec{i}sin(x-vt)[/math]
Es necesario redefinir el parámetro [math]\sigma[/math]:
[math]\sigma[/math] = [math]\lambda\bigtriangledown \cdot \vec{u}Id+2\mu\epsilon[/math] = [math]\lambda\cdot \frac{2}{5}cos(x-vt)\cdot Id+2\mu\frac{2}{5}cos(x-vt)[/math] = [math](\lambda+2\mu)\frac{2}{5}cos(x-vt)[/math]
Procedemos a calcular la divergencia de [math]\sigma[/math]
[math]\bigtriangledown \cdot \sigma = \begin{pmatrix} (\lambda+2\mu)\frac{2}{5}cos(x-vt) & 0 & 0 \\ 0 &(\lambda+2\mu)\frac{2}{5}cos(x-vt) & 0 \\ 0 & 0 & (\lambda+2\mu)\frac{2}{5}cos(x-vt)\end{pmatrix}[/math] = [math]\begin{pmatrix}-(\lambda+2\mu)\frac{2}{5}sen(x-vt)\\ 0\\ 0\end{pmatrix}[/math]
Ahora realizamos [math]\frac{\partial ^2\vec{u}}{\partial t^2}[/math]
[math]\frac{\partial ^2\vec{u}}{\partial t^2}[/math] = [math](-\frac{2}{5}v^2sen(x-vt))\vec{i}[/math]
Entonces:
[math]\vec{F}[/math] = [math]\frac{\partial ^2\vec{u}}{\partial t^2}[/math] - [math]\bigtriangledown \cdot \sigma[/math] = [math](-\frac{2}{5}v^2sen(x-vt))\vec{i}[/math] + [math]((\lambda+2\mu)\frac{2}{5}sen(x-vt))\vec{i}[/math]
Tras aplicar la condición de que [math]\vec F=0[/math] se procede a despejar el parámetro deseado, es decir, [math]v[/math]. Por lo tanto:
[math]v^2[/math] = [math](\lambda+2\mu)\vec{i}[/math]
v= [math]\sqrt(\lambda+2\mu)\vec{i}[/math]
11.1 Cálculo de la velocidad de propagación
Si la onda es transversal,[math]\:\vec{a}=\frac{2}{5}\vec{j}\:[/math], la velocidad de propagación se define como:
Sabiendo que [math]\sigma[/math]: no varia, procedemos al cálculo de [math]\frac{\partial ^2\vec{u}}{\partial t^2}[/math] Entonces: [math]\vec{F}[/math] = [math]\frac{\partial ^2\vec{u}}{\partial t^2}[/math] - [math]\bigtriangledown \cdot \sigma[/math] = [math](-\frac{2}{5}v^2sen(x-vt))\vec{j}[/math] + [math]((\lambda+2\mu)\frac{2}{5}sen(x-vt))\vec{i}[/math]
Aplicamos la condición [math]\vec F=0[/math] y procedemos a calcular [math]v[/math]. Por lo tanto: [math]v^{2}\vec{j}=(\lambda+2\mu )\vec{i} [/math]
La única solución que verifica la igualdad vectorial es [math]\vec{0}[/math]
Las ondas longitudinales y transversales no viajan a la misma velocidad por el mismo medio, tal y como queda demostrado.
12 Módulo de desplazamiento vertical en el punto (1/2,1)
La onda que se trata en este trabajo es una onda longitudinal, es decir, que la onda se desplaza en el mismo sentido qeu la amplitud, por ende, se puede decir con certeza que no se obtendrán desplazamiento en otro eje disntinto al eje en el que se transmite dicha onda. No obstante, se procede a la demostración:
Adoptando el campo de desplazamientos: [math]\vec u (x,y,t)=\frac{2}{5}\vec i sen(x-vt)[/math]
Fijando el punto (x,y)=(1/2,1) se tiene: [math]\vec u (1/2,1,t)=\frac{2}{5}\vec i sen(1/2-vt)[/math] donde el parámetro [math]t[/math], tiempo, pertence al intervalo [math][0,10][/math]
En la dirección de [math]\vec j [/math]se obtiene que: [math]\vec u \cdot \vec j = \frac{2}{5}sen(1/2-vt)\vec i \cdot \vec j = 0[/math] como se esperaba por tratarse de una onda longitudinal.
En la dirección de [math]\vec i [/math]se obtiene que: [math]\vec u \cdot \vec j = \frac{2}{5}sen(1/2-vt)\vec i \cdot \vec i = \frac{2}{5}sen(1/2-vt)[/math]
Es decir, como conclusión se obtiene lo esperado, al ser una onda longitudinal en la que su dirección de propagación es el mismo que la amplitud, se producen desplazamientos horizontales y no verticales.








